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ABSTRACT 


A FORTRAN IV computer program was written to obtain a local detailed solution 
around a leading or trailing edge or in a slot region for compressible, subsonic, non- 
viscous flow on a blade-to-blade surface between turbomachine blades. This program 
allows a coarse- mesh solution for an entire blade-to-blade region to be magnified in a 
small rectangular region. The results include detailed surface velocities, velocity 
magnitude and direction, and stream-function values throughout the magnified region. 
The method is based on the stream function and uses the iterative solution of nonlinear 
finite-difference equations. 
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FORTRAN PROGRAM FOR CALCULATING VELOCITIES IN A MAGNIFIED REGION 
ON A BLADE-TO-BLADE STREAM SURFACE OF A TURBOMACHINE 
by Theodore Katsanis and William D. McNally 
Lewis Research Center 

SUMMARY 

A FORTRAN IV computer program was written to obtain a local detailed solution 
around a leading or trailing edge or in a slot region for compressible, subsonic, non- 
viscous flow on a blade -to -blade surface between turbomachine blades. This program 
allows a coarse-mesh solution for an entire blade -to -blade region to be magnified by a 
chosen magnification factor in a small rectangular region. 

The program input requires information obtained from a less detailed solution from 
one of three other FORTRAN programs. These programs have been presented in NASA 
Technical Notes. The output includes detailed surface velocities, velocity magnitude 
and direction, and stream-function values throughout the magnified region. 

The method is based on the stream function with the solution of the simultaneous, 
nonlinear, finite-difference equations being obtained by two major levels of iteration. 

The inner iteration consists of the solution of simultaneous linear equations by successive 
over relaxation, using an estimated optimum overrelaxation factor. The outer iteration 
then changes the coefficient of the simultaneous equations to compensate for compres- 
sibility. 

This report includes the FORTRAN IV computer program with an explanation of the 
equations involved, the method of solution, and the calculation of velocities. Numerical 
examples have been included to illustrate the use of the program and to show the results 
which are obtained. 


INTRODUCTION 

In the design of blade rows for turbines and compressors , it is desirable to obtain 
the velocity distribution through the passage and particularly over the blade surfaces. 
The authors have published computer programs (refs. 1 to 3) for obtaining this type of 



solution for single and tandem blade rows . 

With these programs, however, it is not always possible to obtain sufficient detail 
on some critical parts of the blade surfaces. These programs give an approximate solu- 
tion for velocities only at the mesh points of a finite-difference grid. Due to storage 
requirements on the computer, grid spacing may be too large to give the desired detail 
around small leading- or trailing- edge radii or within slot regions. And it is in these 
regions where geometry, and thus velocities, change most rapidly. 

For these reasons a computer program has been written to obtain a solution on a 
fine mesh in a small part of the blade-to-blade region. The method used is similar to 
that used by Kramer (ref. 4). A small rectangular region of the solution obtained by 
either 2DCP (ref. 1), TURBLE (ref. 2), or TANDEM (ref. 3) can be magnified by a cho- 
sen factor using MAGNFY, the program described herein. The input and output are sim- 
ilar to 2DCP, TURBLE, and TANDEM, with additional input required. The additional 
input is obtained from the output of 2DCP, TURBLE, or TANDEM. 

This report includes the FORTRAN IV computer program that was developed with 
explanationof the input required. An axial-flow turbine rotor slot and the tip of a mixed- 
flow impeller have been analyzed to illustrate the use of the program. 

This report is organized so that the engineer desiring to use this program needs to 
read only the sections MATHEMATICAL ANALYSIS, NUMERICAL EXAMPLES, and 
DESCRIPTION OF INPUT AND OUTPUT. The necessary information of interest to a pro- 
grammer is contained in the sections DESCRIPTION OF INPUT AND OUTPUT and PRO- 
GRAM PROCEDURE. 

A MAGNFY source deck on tape is available from COSMIC (Computer Software 
Management and Information Center), Computer Center, University of Georgia, Athens, 
Georgia 30601. The program number is COSMIC number LEW-10789. 


SYMBOLS 


A 

b 


m 

R 


coefficient matrix, (eq. (A7), ref. 1) 

stream- channel thickness normal to meridional streamline, meters 
/k., 


constant vector, 


, (eq. (A7), ref. 1) 


\kJ 

meridional streamline distance, meters 
gas constant, joule/(kg)(°K) 
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r radius from axis of rotation to meridional stream -channel mean line, meters 

T temperature, °K 

u stream function 

W fluid velocity relative to blade, meters/sec 

w mass flow per blade flowing through stream channel, kg/sec 

/3 angle between relative velocity vector and meridional plane, rad 

y specific-heat ^ratio 

6 relative angular coordinate, rad 

X prerotation (rVaj , meters^/sec 

' w/ in 

Q 

p density, kg/meters 

O overrelaxation factor, (eq. (A8), ref. 1) 

co rotational speed, rad/sec 

Subscripts: 

in inlet or upstream 

le leading edge 

m component in direction of meridional streamline 

te trailing edge 

6 tangential component 

Superscript: 

’ absolute stagnation condition 


MATHEMATICAL ANALYSIS 

It is desired to determine the flow distribution over the leading or trailing edge of a 
turbomachine blade or through the slot of a tandem or slotted blade. The stream func- 
tion is used for the analysis. The basic assumptions and equations are given in refer- 
ences 1 and 3. The only difference in this analysis from that of references 1 and 3 is in 
the boundary conditions. For the MAGNFY program, the value of the stream function 
must be given for the entire boundary of the region considered. These values are deter- 
mined on a coarse mesh by 2DCP, TURBLE, or TANDEM. MAGNFY then interpolates 
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these values to obtain boundary values of stream functions on a finer mesh. These 
boundary conditions determine a solution to the stream function (eq. (1), ref. 1). The 
numerical solution is determined by using finite-difference equations, as described in 
appendix A of references 1 and 3. 


NUMERICAL EXAMPLES 

Two numerical examples are given to illustrate the use of the program and to show 
the type of results which can be obtained. The first example is the slot region of a tan- 
dem axial-flow turbine rotor blade, and the other is the trailing edge of a mixed-flow 
impeller. 


Leading Edge of Rear Blade of Tandem Blade Turbine Rotor 

Flow about the leading edge of the rear blade of a tandem axial-flow turbine rotor 
blade (ref. 5) has been analyzed to illustrate the use of MAGNFY. The entire blade was 
first analyzed by using TANDEM with the mesh size shown in figure 1 . Due to computer 
storage limitations, this was as fine a mesh as could be obtained with TANDEM. How- 
ever, more detail was desired for velocities between adjacent mesh points on the leading 
edge of the tandem blade. Therefore, MAGNFY was used in order to reduce the mesh 
spacing in the region around the leading edge by a factor of 8, as shown in figure 2. 

The input for this case is given in table I. This includes most of the input necessary 



Figure 2. - Magnified region of figure 1 with reduced 
mesh size. 



Figure 1. - Tandem axial turbine with region to be magnified. 
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for TANDEM, in addition to the stream-function values on the vertical and horizontal 
boundaries of the magnified region. These stream-function values and their coordinates 
were obtained as output from the TANDEM program. 

Blade-surface velocities from both TANDEM and MAGNFY are plotted in figure 3. 
The original velocities obtained by using the coarse mesh are denoted by circles in fig- 
ure 3. The MAGNFY output is plotted as a solid line. As shown in figure 3, the velocity 
peak on the suction surface is much higher than indicated by the coarse-mesh solution. 
This illustrates the need for a finer grid in some parts of the solution region since the 
velocities denoted by circles do not define the velocity adequately. 

The execution time for this example was 10 minutes on the direct -coupled IBM2- 
7094-7 044 computer . 


O TANDEM program 



Figure 3. - Comparison of velocities from 
coarse-mesh and fine-mesh solutions. 


Trailing Edge of Mixed-Flow Impeller 

Flow about the trailing edge of the main blade of a mixed-flow impeller with splitter 
blades (ref. 6) has also been analyzed using MAGNFY. The entire impeller was origin- 
ally analyzed with TANDEM (ref. 3). The impeller profile (in the meridional plane) is 
shown in figure 4(a). Figure 4(b) shows the blade-to-blade region with the coarse mesh 
used in the TANDEM run. The region to be magnified about the trailing edge of the main 
blade is indicated by heavy lines in figure 4(b). MAGNFY was used to reduce the mesh 
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TABLE I. - MAGNFY INPUT FOR LEADING EDGE OF REAR BLADE OF TANDEM BLADE TURBINE ROTOR 


SHIFTED TANOEM AXIAL TURBINE ROTUR - SMALL 
NOBL 

SLOT AREA 





c 

GAM 

1.4000000 
BETA I 
48.C00000 
MBI MBO MB 12 
10 32 24 

AR 

287.05300 
BET AO 

-48.000000 
M602 MM NBBI 
49 58 20 

TIP 

288. 15000 
CHOROF 

0 • 2847000E-01 
NBL NRSP 
76 2 

RHQIP 

1.2250000 

STGRF 

0. 21 33300E-01 

WTFL 

0.2715000 

CHDROR 

0.3183000E— 01 

WTFLSP 

0.6970000E— 01 
STGRR 

-0.5647000E-01 

OMEGA 

-0 

ML ER 

0.17730 00 E-Ol 

ORF 

1.7300000 

THLER 

0.2820000E-02 


BLADE SURFACE 1 — UPPER SURFACE - FRONT BLADE 

RI1 ROl BETH BETOI SPLNOL 

0.7620000E-03 0. 38 10000E-03 50.000000 -29.400000 7.0000000 

MSP L ARRAY 

-0 0. 2570000E— 02 0 . /650000E-02 0. 1527000E-0! 0.2035000E-01 0 . 2543000E-01 -0 

THSP1 ARRAY 

-0 0. 92500 00E- 02 0 . 21 1 8000E-01 0 . 29B8000E-01 0.3020000E-01 0. 2643000E-01 -0 


BLADE SURFACE 2 — LOWER SURFACE - FRONT BLADE 

R 1 2 R02 BET 1 2 BET02 SPLN02 

0.7620000E-03 0 . 38 10000E-03 25.000000 -6.9000000 5.0000000 

MSP2 ARRAY 

-0 0. 7650000E— 02 0.2035000E-01 0 . 2543000E-01 -0 

THSP2 ARRAY 

-0 0.7140000E-02 0 . 2039000E-01 0 . 2094000E-01 -0 


BLADE SURFACE 3 -- UPPER SURFACE - REAR BLADE 

R 1 3 R03 BET 1 3 BET03 SPLN03 

0. 7620000E-03 0 . 38 10000E-03 13.000000 -48.800000 5.0000000 

MSP J ARRAY 


-0 

0. 5160000E-02 

0, 1278000E-01 

0.2294000E-01 

-0 

THSP3 

-0 

ARRAY 

0.408 00 OOE- 02 

-0.2400000E-03 

-0.2651000E— 01 

-0 


BLADE SURFACE 

4 — LOWER SURFACE - REAR BLADE 


R 1 4 

R04 

BET 1 4 

BET04 

SPLN04 

0. 7620000E-0 3 

0 . 38 10000E— 03 

-4.2000000 

-42.500000 

5.0000000 

MSP4 ARRAY 

-0 

0.5 160000E-02 

0. 1278000E-01 

0.2040000E-01 

-0 

THSP4 ARRAY 

-0 

-0.439 00 OOE— 02 

-0. 1388000E-01 

-0 • 293 3 000 E— 01 

-0 


MR ARRAY 

-l.OCOOOOO 1.0000000 

RMSP ARRAY 

0. 3238 300 0.3238500 

BESP ARRAY 

0. 1143000 0.1143000 

BLOAT AANOK ERSOR STRFN INTVL SURVL 

1 1 2 2 2 3 



MBDYF MBDYL ITF I TL MAGFAC 

22 27 -1 A 8 

LAMBDA 
27.587140 


KBDRY NSP 
1 6 
BVIN ARRAY 




0. 1519700E-01 
UBVIN ARRAY 

0. 1646400E-01 

0.1773000E-01 

0.1907200E-01 

0.2562700 

-0.2741800 

-0.2883600 

-0.2929900 


KBDRY NSP 
2 6 
BVIN ARRAY 

0. 1519 700E-0 1 0.1646400E-01 0. 1773000E-01 0. 1907200E-01 

UBV IN ARRAY 

-0. 1300000F-03 -0.1841000E-01 -0 . 3630000E-01 -0.5292000E-01 


KBDRY NSP 
3 6 

BVIN ARRAY 
-0.4133700E-02 0 

UBV I N ARRAY 

-0.2562700 -0.2123600 


KBDRY NSP 
4 2 

BVIN ARRAY 

— 0.41337 00 E-0 2 -0.8230000E-03 

UBVIN ARRAY 

-0.2874300 -0.2567200 


KBDRY NSP 
4 4 

BVIN ARRAY 




0.6708000E-02 
UBVIN ARRAY 

0.8267400E-02 

0. 1240100E-01 

0. 1653500E-01 

-0.2567200 

-0.2261100 

-0.1478300 

-0. 7366000E-01 


0.4133700E-02 0 . 8267400E-02 

-0.1671100 -0.1146800 


-7 


0.2041 500E-01 0.2175700E-01 

-0.292D900 -0.2874300 

0.2041 500E-01 0. 2175700E-01 

-0.6547 00 0E-01 -0 . 7366000E-01 

0.1240100E-01 0. 1653500E-01 

-0.5821 OOOE-Ol -0. 1300000E-03 





Angular coordinate, 6, radians 



Meridional streamline distance, m, meters 
(b) Blade-to-blade surface, showing coarse grid used in TANDEM program. 
Figure 4. - Concluded. 
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TABLE H. - MAGNFY INPUT FOR TRAILING EDGE OF MIXED- FLOW IMPELLER 


MIXED FLOW IMPELLER (NASA TN D-1186) - TRAILING EDGE, MAIN BLADE 
NOBL 


2 



GAM 

AR 

TIP 

RHOIP 

WTFL 

WTFLSP 

1 . 

5C00000 

1000.0000 

1000000.0 

1.0000000 

0.3042000E-02 

0. 1351600E-02 


BET A I 

BETAO 

CHORDF 

STGRF 

CHORDR 

STGRR 

-84 

.080000 

-43.000000 

0. 1055500 

-2.6290000 

0.5664000E— 01 

-0.6649000 

MBI 

MBO MBI2 

MB02 MM NBBI 

NBL NRSP 




LO 

47 28 

47 57 28 

8 18 





OMEGA 

796,00000 

MLER 

0.48910C0E-01 


ORF 

1.8440000 

THLER 

-2.3434000 


BLADE SURFACE 1 — UPPER SURFACE - FRONT BLADE 


R 1 1 

R01 

BETH 

BET01 

SPLN01 


0.9140000E-03 
MSP l ARRAY 

0. 18460006-02 

-80.000000 

-49.000000 

6.0000000 


0 

THSP1 ARRAY 

0.1214000E-Q1 

0.2651000E-01 

0.4766000E-01 

0.7360000E-01 

0 

0 

-0.6250000 

-1.2330000 

-1.8182000 

-2.2750000 

0 


BLADE SURFACE 

2 — LOWER SURFACE - FRONT BLADE 



R I 2 

R02 

BET I 2 

BET02 

SPLND2 


0.9 140000E-03 

0 • 1 84600 OE- 02 

-83.000000 

-41.500000 

6.0000000 


MSP2 ARRAY 






0 

0.788 00 OOE- 02 

0.2004000E-01 

0.40060 00E-01 

0.6828000E-01 

0 

THSP2 ARRAY 






0 

-0.6310000 

-1.2310000 

-1.8206000 

-2.2954000 

0 

BLADE SURFACE 

3 — UPPER SURFACE - REAR BLADE 



R 1 3 

R03 

BET 1 3 

BET03 

SPLN03 


0. 1328000E-02 

0. 1753000E-02 

-60.500000 

-51.500000 

6.0000000 


MSP3 ARRAY 






0 

0. 1307000E-01 

0.2552000E-01 

0.4172000E— 01 

0.5280000E-01 

0 

THSP3 ARRAY 






0 

-0. 1670000 

-0.3370000 

-0.5262000 

-0.6269000 

0 


BLADE SURFACE 4 — LOWER SURFACE - REAR BLADE 


R 1 4 

R04 

BETI4 

BET04 

SPLN04 

0. 1328000E-02 

0. 1753000E-02 

-63.000000 

-40.500000 

5.0000000 

MSP4 ARRAY 

0 

0. 1073000E— 01 

0 .249 30 OOE- 01 

0.4172000E-01 

0 

THSP4 ARRAY 

0 

•0.2010000 

-0.4070000 

-0.5818000 

0 


MR ARRAY 

-0.3124000E-01 

-0.151 40 00E-01 

0.2500000E-03 

0.1065000E-01 

0. 1853000E-01 

0. 265 1000E— 01 

0.34600 C0E-01 

0.4281000E-01 

0.51 15000E-0 1 
0. 1272600 

0 • 5964000E-01 
0.1407300 

0.6828000E-01 

0.7709000E-01 

0.8607000E— 01 

0 • 9524000E-01 

0.1046100 

0.1141700 

RMSP ARRAY 

0. 7586000E-0 1 

0. 7662000E-01 

0.7874000E-01 

0.8091000E— 01 

0.8294000E-01 

0.8531000E-0L 

0.8802000E— 01 

0.9108 00 OE -01 

0.9447000E-01 

0.1360200 

0 • 9820000E-01 
0.1448700 

0.1022800 

0.1067400 

0.1114600 

0.1165600 

0.1220000 

0.1277800 

BESP ARRAY 

0. 1053300E-0 1 

0. 1004500E— 01 

0 . 8724000E-02 

0 • 7420000E-02 

0.6316000E-02 

0 • 5354000E-02 

0.45 320 OOE— 02 

0. 3831000E-02 

0. 3235D00E-02 
0.8250000E-03 

0.2728000E-02 

0.7240000E-03 

0.2299000E-02 

0. 1936000E-02 

0. 1629000E-02 

0. 1370000E-02 

0. 1151000E-02 

0.9790000E-03 



BLDAT AANDK ERSOR STRFN 
12 2 2 


MBOYF HBOYL ITF ITL 
44 50 -97 -90 

LAMBDA 

-0.9835000E-02 


KBDRY NSP 
1 7 

BVIN ARRAY 

0. 9660 700E-0 1 0 . 9958800E-01 

0.1025700 

0.1055500 

0.1085300 

0.1115100 

0.11449C0 

UBVIN ARRAY 

-0.2089400 -0.1757100 

-0.1440400 

-0.1137000 

-0.8342000E— 01 

-0.5251000 £-01 

-0.2076000E-01 

KBDRY NSP 
2 7 

BVIN ARRAY 

0.9660700F-01 0. 9958800E-01 

0.1025700 

0.1055500 

0.1085300 

0.1115100 

0.1144900 

UBV1N ARRAY 

0.1517000E-01 0.4847000E-01 

0.8207000E-01 

0.1159000 

0.1498600 

0.1839400 

0.2181900 

KBDRY NSP 
3 6 

BVIN ARRAY 

-2.7208400 -2. 6927900 

-2.6647400 

-2.6366900 

-2.6086400 

-2.5949900 


U8VIN ARRAY 

-0.2089400 -0.1676100 

-0.1226400 

-0 . 7412000E-01 

-0. 2 446 00 OE -01 

0 


KBDRY NSP 

3 2 
BVIN ARRAY 

-2.5380100 -2.5244900 

UBVIN ARRAY 

0 0. 151 7000E-01 

KBDRY NSP 

4 8 
BVIN ARRAY 

-2.7208400 -2.6927900 

-2.6647400 

-2.6366900 

-2.6086400 

-2.5805900 

-2.5525400 

UBVIN ARRAY 

— 0.2077000E-0 1 0. 1455000E-01 

0.4942000E-01 

0.8343000E-01 

0.1166400 

0.1499700 

0.1837800 


INTVL SURVL 
2 3 


MAGFAC 

5 


-2.5244900 

0.2181900 



Distance along surface, cm 

Figure 6. - Comparison of velocities near 
trailing edge of mixed-flow impeller. 


size in this region by a factor of 5. The reduced mesh is shown in figure 5. 

The input for this example is given in table II . It includes the original TANDEM 
input plus stream-function boundary values about the magnified region. These were 
obtained from the output of TANDEM. The process by which this input was obtained is 
explained later in the section Example of Preparing Input. 

The blade-surface velocities from TANDEM and MAGNFY are plotted in figure 6. 
The velocities obtain by TANDEM (coarse mesh) are not accurate near the trailing-edge 
radius, and the magnitude of the peak is not shown accurately. 

The execution time for this example was 2 minutes on the direct -coupled IBM2- 
7094-7044 computer. 
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DESCRIPTION OF INPUT AND OUTPUT 


The computer program requires as input the same input as was used in either 2DCP 
(ref. 1), TURBLE (ref. 2), or TANDEM (ref. 3) plus the stream-function values along 
the boundary of the region to be magnified. These stream-function values are generally 
obtained from the output of either 2DCP, TURBLE, or TANDEM (refs. 1 to 3). 

Output obtained from MAGNFY includes velocity magnitude and direction at all 
interior mesh points in the region, blade-surface velocities, and stream-function values 
throughout the region. 


Instructions for Preparing Input 

The first step in obtaining input for MAGNFY is to obtain the usual coarse-mesh 
solution from either 2DCP (ref. 1), TURBLE (ref. 2), or TANDEM (ref. 3). If TURBLE 
or TANDEM are used, their input forms part of the input for MAGNFY. If 2DCP is used, 
modifications must be made to its input to make it appropriate for MAGNFY. These 
modifications make the 2DCP input agree with input for TURBLE, and are explained in 
the following section. 

The remainder of the input for MAGNFY (beyond the input for 2DCP, TURBLE, or 
TANDEM) consists of coordinates and stream-function values obtained from the output 
of one of these three programs. Figure 7 shows all the input for MAGNFY beyond that 
required for 2DCP, TURBLE, or TANDEM. 

The boundary value input for MAGNFY could be determined from some method other 
than either 2DCP, TURBLE, or TANDEM. In this case the input which would ordinarily 
have been used with these programs must be determined as explained in references 1 or 2 
(single blade) or reference 3 (tandem or slotted blade). 


Modification of 2DCP Input 

If the user desires to magnify a solution obtained with 2DCP, the 2DCP input must 
be rearranged as if it were to be run on TURBLE before it can be used with MAGNFY. 
Some of the 2DCP input variables have the same names as MAGNFY variables, but some 
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MSP2 ARRAY 

w/////#m 

THSP2 ARRAY 

y////////sWA 

MR ARRAY 

I 

I 

RMSP ARRAY 


BESP ARRAY 



Figure 8. - Input form for TURBLE (ref. 2). Card column numbers appear at top. 

do not. Table HI lists the 2DCP and MAGNFY variables which have the same meaning 
but different names. Figure 8 and table III show the user how to rearrange his 2DCP 
input to make it compatible with MAGNFY. The first card of 2DCP input (the GAM card) 
must be modified for MAGNFY by shifting OMEGA and W 10 spaces to the right (see 
fig. 8). The second card should have the inlet and outlet flow angles (BETAI and BETAO) 
placed before CHORD and STGR. Also BETAI and BETAO have been redefined to be the 
flow angles at the leading and trailing edges , instead of at upstream and downstream 
boundaries. The third card contains information obtained from the fourth card for 2DCP. 
Once again, the position of variables on the card and the relation between 2DCP and 
MAGNFY variables can be seen from figure 8 and table III. The information on the third 
2DCP card must be placed on two cards for MAGNFY. These cards are placed directly 
above the two sets of m- and 0-coordinates (MSP1, 2 and THSP1, 2) for the two blade 
surfaces. Each of these cards for MAGNFY contains inlet and outlet radii, tangency 
angles, and number of spline points for one of the two blade surfaces. Finally, the cards 
containing m- and 0-coordinates are unchanged between 2DCP and MAGNFY. 


To 
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TABLE m. - 2DCP AND MAGNFY VARIABLES WITH 
SAME MEANING BUT DIFFERENT NAMES 


2DCP variable 


Corresponding MAGNFY variable 


w 

ORF 

RI 

RI1 and RI2 

RO 

ROl and R02 

ALUI 

BETIl 

ALLI 

BETI2 

ALUO 

BETOl 

ALLO 

BET02 

MXBI 

MBI 

MXBO 

MBO 

i 

MX 

MM 

NUSP 

SPLNOl (real) 

NLSP 

SPLN02 (real) 

MU 

MSP1 

XSPU 

THSP1 

ML 

MSP2 

XSPL 

THSP2 


Choosing Magnified Region 

The region where magnification is desired is usually located around a leading or 
trailing edge or about a tandem blade slot. Therefore, the region generally includes 
portions of both lower and upper blade surfaces. This is indicated by the heavily 
outlined region of figure 9 and by the similar region for the mixed-flow impeller example 
(fig. 4(b)). However, input must be given to MAGNFY as though the region was entirely 
located about the lower blade. (See the dashed portion of the region in fig. 9. ) This con- 
dition, of course, results in a magnified region which is partially outside of the original 
2DCP, TURBLE, or TANDEM region. In the case of the leading edge of the rear blade 
of a tandem blade, the magnified region may lie completely outside of the original 
TANDEM region. The region may contain at most one leading and one trailing 
edge. 

The fact that the magnified region is restricted to the lower blade implies that, once 
it is drawn about the lower blade, no part of it may include any of the upper blade. 
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Input 


Only the additional input beyond that required for 2DCP, TURBLE, or TANDEM will 
be described in detail herein. All integers are in a five-column field and must be right 
adjusted. All numbers with a ten-column field are real numbers and must have a decimal 
point punched on the data card. 

The first input data card (fig. 7) is a label card containing any desired identification 
label. The next card has either NOBL=l (single blade) or NOBL=2 (tandem or slotted 
blade) in column 5. This is followed by the input data from either TURBLE or TANDEM, 
or modified input from 2DCP. This input consists of all the input cards from the first 
one starting with GAM up to the last geometry input card containing values of BESP. 

This input remains unchanged, except for the overrelaxation factor (ORF), which should 
be recalculated for the MAGNFY program; that is, it should be set equal to zero again 
for the initial MAGNFY run on a set of data. 

The next input card has variables (BLDAT to SURVL) used to indicate what output 
is desired. These variables are used in the same way as in TURBLE or TANDEM, 
except for the ommission of SLCRD, which is not required in MAGNFY. 
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The remaining input variables (see fig. 7) are as follows: 

MBDYF index IM of vertical mesh line which is to be left boundary of magnified 
region (see fig. 9) 

MBDYL index IM of vertical mesh line which is to be right boundary of magnified 
region 

ITF index IT of horizontal mesh line which is to be lower boundary of magnified 

region (See previous section for explanation of how to choose magnified 
region. ) 

ITL index IT of horizontal mesh line which is to be upper boundary of magnified 

region 

MAGFAC magnification factor (If MAGFAC = n, one mesh square of original coarse 

2 

mesh will contain n squares of smaller mesh. ) 

LAMBDA value of prerotation X at inlet (LAMBDA is given as part of the output for 
2DCP, TURBLE, or TANDEM. ) 

KBDRY indicates which boundary is referred to on input cards which follow it: 
KBDRY = 1 , low e r bo undary 
KBDRY=2, upper boundary 
KBDRY=3, left boundary 
KBDRY=4, right boundary 

NSP number of stream -function values given for a particular boundary on BVTN 

and UBVIN cards which follow it 

BVTN boundary coordinates (m or 0 ) for segment of boundary indicated by KBDRY 

(These coordinates should correspond to original coarse-mesh lines, 
except for possibly the first and last points, which could fall on a blade 
surface. ) 

UBVIN stream-function values corresponding to BVTN 

The variables from KBDRY to UBVTN are given for each segment of the boundary. After 

all boundary values are given for each segment of the four boundaries, one blank card 

(or a card with zeros for KBDRY and NSP) is added to signal the end of the input data for 

a particular case. 

Example of Preparing Input 

The second numerical example of this report (p. 5) dealt with solving for detailed 
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velocities about the trailing edge of the main blade of a mixed-flow impeller with splitter 
vanes (ref. 6). This section illustrates, in detail, how the additional MAGNFY input 
(beyond the normal TANDEM input) for that example was obtained. This input is given 
in table II. 

MAGNFY is intended to be used by those who have run 2DCP, TURBLE, or TANDEM, 
and who desire a more detailed solution about some critical region on the blade. For 
the impeller example, that region is the trailing edge of the main blade as shown in fig- 
ure 10. Notice that most of this region lies about blade surface 1, but that a portion of 
it is located at the end of blade surface 2. In most cases, the region to be analyzed will 
be divided in this way. However, input must be given to MAGNFY as though the region is 
entirely located about the lower blade . The way this is done is illustrated in the follow- 
ing paragraph. 

The user should draw a magnified picture of the region for which a detailed solution 
is desired (fig. 11). This rectangular picture should extend three or four mesh lines 
(coarse mesh) in all directions from the point at which most detail is desired. The 
coarse mesh should be numbered with IM and IT grid line values. The ITV array of 
TURBLE or TANDEM can be used in the drawing of this sketch. The boundaries parallel 
to the IM axis are defined as boundaries 1 and 2 of this sketch, and those parallel to the 
IT axis as 3 and 4, as indicated in figure 11. 

At this point, some of the input to MAGNFY can be obtained. The first and last 
values on the IM and IT axes are called MBDYF, MBDYL, ITF, and ITL. For the impel- 
ler example, these values are 44, 50, -97, and -90, respectively. 


Boundary 2 
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Boundary 1 

Figure 11. - Magnified region for input example. 



MAGFAC can also be chosen at this time. A number between 5 and 8 is typical, and 
in this case 5 was used. MAGFAC must be chosen so that the resulting grid has less 
than 2000 mesh points in it. 

LAMBDA is the next required input. It is obtained directly from the output of either 
2DCP, TURBLE, or TANDEM. 

The remainder of the MAGNFY input consists of geometrical (BVIN) and stream- 
function (UBVIN) boundary values for the coarse -mesh boundary points on the four bound- 
aries of the magnified region (fig. 11). The boundaries 1 to 4 should be entered in order, 
giving values from left to right on boundaries 1 and 2 and from bottom to top on boundaries 
3 and 4. 

The blade surfaces always intersect some boundaries of the region. When this 
occurs, the resulting sections of boundaries should be entered separately. In this exam- 
ple, boundary 3 is divided into two parts. The first has six points (IT = -97 to IT - -93 
plus the blade-surface point (point a), fig. 11). The second has two points (point b and 
IT = -90). 

For each section of a boundary, four items of input are needed: KBDRY, NSP, the 
BVIN array, and the UBVIN array. KBDRY identifies by number the boundary for which 
data are given, and NSP is the number of points given on a section of that boundary. For 
the impeller example, NSP = 7 on boundaries 1 and 2; NSP = 6 and NSP = 2 on boundary 3; 
and NSP = 8 on boundary 4. 

For boundaries 1 and 2, BVIN is obtained from the Stream Sheet Coordinates and 
Thickness Table for both TURBLE and TANDEM. A portion of the table for this example 
is reproduced in table IV. The meridional coordinates for BVIN from IM = 44 to IM = 50 
are circled. 


TABLE TV. - STREAM SHEET COORDINATES AND THICKNESS TABLE 


IM 

M 

R 

SAL 

B 

0B/DM 

I 

-0.24455E-01 

0. 76028E-01 

0.3899OE-0L 

0. 10422E-01 

— 0.251 97 E— 01 

2 

-0.21738E-01 

0. 76152E-01 

0.52336E-01 

0. 10342E-01 

-0.33504E-01 

3 

-0. 19021E-01 

0.76314E-01 

0.67I78E-01 

0.10239E-01 

-0.42743E-01 

4 

-0.16303E-0L 

0.765I9E-01 

0.83516E-OI 

0. 10109E-01 

-0.52913E-01 

5 

-0. I3586E-01 

0. 76769E-01 

0.10100 

0. 99507E— 02 

— 0.63772E— 01 


39 

0.81 702 E— 01 

0.10914 

0.52474 

0. 17711E-02 

— 0.34060 E— 01 

40 

0. 84683E— 01 

0.11071 

0.53430 

0.16727E-02 

-0.31962E-01 

41 

0 . 8 7664E— 0 1 

0.11233 

0.54631 

0. 15804E— 02 

-0.29980E— 01 

42 

0 • 90645E— 0 1 

0.11397 

0.55690 

0. 14938E-02 

-0.28176E-01 

43 

0 • 93626E— 01 

0.11564 

0.56553 

0. 14122E-02 

-0. 26560E— 0 1 

44 

0 • 96607 E— 01 

0. 11734 

0.57246 

0. 13352E— 02 

-0.25103E— 01 

45 

0 • 99588E— 0 1 

0.11906 

0.57954 

0. 12626E-02 

— 0.23590 E— 01 

46 

0.10257 

0.12080 

0.58712 

0.11947E-02 

-0.21979E— 01 

47 

0.10555 ■ 

0.12256 

0.59514 

0. 11317E-02 

— 0.2 0275 E— 01 

48 

0.10853 

0.12434 

0.60279 

0. 10739E-02 

-0. 18520E-01 

49 

0.11151 

0.12615 

0.60980 

0. 10213E-02 

-0. 16728E-01 

50 

0.11449 

0.12798 

0.61617 

0.97415E— 03 

-0 . 14901 E— 0 1 

51 

0.11747 

0.12983 

0.62238 

0. 93233E— 03 

-0.13193E-01 

52 

0.12046 

0.13169 

0.62878 

0.89527E— 03 

-0.11706E-01 

53 

0.12344 

0.13357 

0.63535 

0. 86232E-03 

-0. 10440 E-01 

54 

0.12642 

0.13548 

0.64210 

0.83281E— 03 

-0.93952E-02 

55 

0.12940 

0.13740 

0.64879 

0. 8061 1 E— 03 

— 0.8 5392 E— 02 
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TABLE V. - THETA COORDINATES OF 
HORIZONTAL MESH LINES 


IT 

THETA 

-107 

-3.00136 

-106 

-2.97329 

-105 

-2.96526 

-106 

-2.91719 

-103 

-2.88916 

-102 

-2. 86109 

-101 

-2.83306 

-100 

-2.80699 

-99 

-2.77696 

-98 

-2.76889 

-97 

-2. 72086 

-96 

-2.69279 

-95 

-2.66676 

-96 

-2.63669 

-93 

-2.60866 

-92 

-2.58059 

-91 

-2.55256 

-90 

-2.52669 

-89 

-2.69666 

-88 

-2.66839 

-87 

-2.66036 

-86 

-2.61229 

-85 

-2.38626 

-86 

-2.35619 


If a blade surface passes through boundary 1 or 2 (it does not in this example), 
the BVTN for the point of intersection of the blade and boundary is obtained from the MH 
array (m- coordinates of intersections of horizontal mesh lines with blade given as output 
from 2DCP, TANDEM, or TURBLE) for the blade surface involved. 

For boundaries 3 and 4, BVTN is obtained from the table of Theta Coordinates of 
Horizontal Mesh Lines for both TURBLE and TANDEM. A portion of this table for this 
example is reproduced in table V. The tangential coordinates for BVIN from IT = -97 
to IT = -90 are circled. 

If a blade surface passes through boundary 3 or 4, the BVIN for the point of inter- 
section of the blade and boundary is obtained from the TV array (0-coordinates of blade 
at vertical mesh lines) for the blade surface involved. In the example, blade surfaces 1 
and 2 pass through boundary 3. A portion of the TV array output from TANDEM for 
surfaces 1 and 2 is given in table VI. The m- coordinate corresponding to IM = 44 is 
circled along with the 0 values, called TV, where the IM = 44 mesh line meets surfaces 
1 and 2. The 0 for blade surface 2 (or 4 on the rear blade of a tandem blade) must 
always have PITCH subtracted from it to bring it down on the same blade as surface 1 
(or 3). 
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TABLE VI. - EXAMPLE OF TV ARRAY OUTPUT FROM TANDEM 
FOR BLADE SURFACES 1 AND 2 

BLADE DATA AT INTERSECTIONS OF VERTICAL MESH LINES WITH BLADES 


BLAOE SURFACE 1 


BLADE SURFACE 2 


M 

TV 

DTDMV 

TV 

DTDMV 

0 

0 

0.10000E 11 

0. 78540 

— 0. 10000E 11 

0.271 72E-02 

-0 . 62034E-0 1 

-69.7931 

0.52933 

-85.1682 

0.54344E-02 

-0.24292 

-63.4793 

0.31771 

-71.3544 

0.81517F-02 

-0.40773 

-57.9565 

0.13743 

-62.0849 

0. 10869F-01 

-0. 55860 

-53.2247 

—0.2 1 730E— 01 

-55.2224 

0.81702E-01 

-2. 37840 

-11.8861 

-1.66869 

-10.7007 

0. 84683E-01 

-2.41299 

-11.3360 

-1.69985 

-10.2016 

0.87664E-0L 

-2.44605 

-10.8536 

-1.72952 

-9.70233 

0.90645F-01 

-2.47777 

-10.4388 

-1.75770 

-9.20294 

0-9 3626E-0 I 

-2.50836 

-10.0916 

-1.78439 

-8.70342 

1 0 • 96607F-0 1 

-2.53801 1 

-9.81206 

1 -1.80959 1 

-8.20376 

0.99588E-0L 

-2.56692 

-9.60010 

-1.83330 

-7.70396 

0. 10257 

-2.59531 

-9.45574 

-1.85548 

-6.36166 

0.10555 

-2.62900 

-0.10000E 11 

-1.84360 

0.10000E 11 


All values of UBVTN are obtained from the table Stream- Function Values. A portion 
of that table for the example has been reproduced in table VII. The boundary values for 
the region of figure 11 have been circled in table VII, and the following paragraph explains 
how they were obtained. 

The table of stream-function values gives u along vertical mesh lines from blade 
to blade. Each mesh line is listed separately, and if a second blade intersects the mesh 
line, the two parts of the mesh line are listed separately. On each part the IT of the 
first mesh point on the line above blade surface 1 or 3 is listed as IT1 in table VII. With 
this information, the proper boundary values for the region can be obtained. 

Values in the shaded portion of the region of figure 11 must be obtained from the 
blade above, as shown in figure 10. 

Along boundary 1 in the example it is desired to obtain u at IT = -97 for IM = 44 
to 50. However, for IM = 44 to IM = 47, u must be obtained from the blade above the 
region in the figure (see fig. 10). Since NBBI (the number of horizontal mesh lines 
between AB and MN) is 28 in this example, u must be obtained at IT = -97 + 28 = -69 for 
these values of IM. 

For IM = 44, ITL = -76 for the upper section of this vertical mesh line. Therefore, 
the eighth value in the row (0. 79106) is the u for IT = -69. To reduce it to correspond 
to IT = -97, the stream-function period (1.0) is subtracted. The input value used is thus 
0. 79106 - 1.0000 = -0. 20894. Likewise, for IM = 45, IT1 = -77; and the ninth value 
(0. 82429) corresponds to IT = -69. Reducing this by 1.0 gives -0. 17571, the value used 
as input. 
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TABLE VII. - STREAM-FUNCTION VALUE TABLE 


IM = 


0. 

0 . 


IN = 


IM = 


0 . 

0 . 


IM = 


E 


0 . 


42 IT1 = 
56370215 0. 
93554489 0. 

43 III = 
01338000 0. 
44875421 0. 

43 I T 1 = 
56540487 0. 
93239887 0. 

44 IT1 = 

OlSlliTVl 0. 

436595451 0. 


IM = 


-74 
59089224 
98355945 
-89 

04592191 
50737296 
75 

59254353 

98098022 

-90 

04746477 

49739579 

-76 

59293610 


0.62004882 0.65130834 0.68482589 0.72076356 0.75924806 0.80028503 0.84364127 0.88884238 


0.08012446 0.11625171 0.15459124 


0.62132733 


0.65190407 


0.68445893 


0.19547091 

0.71922223 


0.23926882 0.28639312 


0.75647087 


0.79649328 


0. 33717524 
0.83945329 


0.39157937 


0.88507573 


44 I\T 1 = 

0.56548788 \ 0. 

10.925880851 \ |0. 97554200 
IM = 45 I T& = -91 

0.01624954 Ip. 048470431 
0.41787726 0.479318821 

IM = 45 IT1 = -77 

0.56370431 0.59194452 

0.91331481 0.96643321 

IM = 46 I T 1 = -92 

0.01692253 0.04918601 

0.39781557 0.45154771 

IM = 46 IT1 = -78 

0.55961629 0.58938974 

0.89759777 0.94536442 

IM = 47 I T 1 = -93 

0.01894281 0.05020627 

0.38282271 0.42824350 


0.08105000 0.11617973 0.15312202 0.19218484 0.23374939 0.27831270 0.32652676 0.37915015 

0.62166080 0.65180898 0.68355786 0.71713170 0.75282977 | 0.791 06467] [o'. 832386291 lo. 87736242 1 

0.36386504 


0.11599385 0.15177523 0.18923291 0.22866775 0.27047658 \ 0.31523652 

0.62101315 0.65109950 0.68237852 0.71504194 0.74932828 0.78557044 


0.08163682 

0.54605063 


10.08207236 

0.51836125) 


0.11586361 0.15076017 0.18695513 

0.61925767 ^0.64969136 0.68091226 0.71309783 

.15012404 0.18537962 0.22180108 0.25954640 


0.22466160 0.26413777 

0.74643017 0.78112448 


0.08260356 

0.47632948 


. 11589979 
0.52582654 





0.86637111 
0.34996708 
0. 855957781 
0.33979683 


1 M 

= f 

111 = “ f V 



/ \ 







0.55809889 

/ 0.58508012 

0.61643707 

0.64747123 

0.67901288 

\ 0.71115961 

0.74405468 

0.77785022 

0.81270609 

0.84879535 


1 0.88629808 

J 0.92530244 

0.96512271 

1.00002798 







IM 

= 4-8 

ITl = -107 










-0.42029791 

-0.38789126 

-0.35575397 

-0. 32350658 

-0.29100066 

-C. 25809307 

-0.22464573 \ 

-0. 19053367 

-0.15565266 

-0.11993713 


1-0.08341560 

] -0.04537519 

-0.00981174 

0.02372140 

0.05218717 

0.08351019 

0.11620351 

0. 149861511 

0.18441523 

0.21990840 


0.25642654 

0.29407255 

0.33295587 

0. 37317996 

0.41480967 

0.45776550 

0.50148795 

0.54395759 



IM 

= 49 

ITl = -107 










-0.394 79662 

-0.36041694 

-0.32688116 

-0.29360276 

-0.26029523 

-0.22677304 

-0. 19289698 

-0. 15856233 

-0.12370409 

-0.08831833 


1-0.05251244 

] -0.01560652 

0.01871358 

0.05240273 

0.08420596 

0.11663655 

0. 14989388 

0.183939301 

0.21877056 

0.25442570 


0.290961571 

0. 32843728 

0.36689877 

0.40635552 

0.44673426 

0.48777 796 

0.52884319 

0.56860353 



IM 

= 50 

ITl = -107 










-0.36652284! 

-0.33123814 

-0.29679292 

-0.2627,4382 

-0.22881144 

-0.19480821 

-0. 16060171 

-0.12610277 

-0.09126864 

-0.05611859 


[-0.02076425 1 10.014548521 

10.049420651 

10.083431531 

|0. 11663859) 

|0. 14997184| 

|0. 183785161 

0.218195711 

0.25326142 

0.28902949 


0.32554135 

0.36282565 

0.40088454 

0.43967149 

0.47905274 

0.51874412 

0.55822743 

0.59672115 



IM 

= 51 

ITl = -107 










-0.33620130 

-0.30048737 

-0.26548242 

-0.23089461 

-0.19649943 

-0.16213188 

-0. 12767380 

-0.09304883 

-0.05822450 

-0.02322172 


0.01187199 

0.04689337 

0.08162532 

0.11590057 

0.14980716 

0.18375354 

0.21801195 

0.25272756 

0.28798253 

0.32382912 


0.36029924 

0.39740183 

0.43511284 

0.47335853 

0.51199066 

0.55075540 

0.58927225 

0.62707399 



IM 

= 52 

ITl = -107 










-0.30426541 

-0.26834252 

-0.23299874 

-0.19804831 

-0.16332335 

-0.12868846 

-0.09404254 

-0.05931953 

-0.02449090 

0.01042976 


0. 0453R360 

0.08027114 

0.11497962 

0.14944728 

0.18374881 

0.21809705 

0.25268143 

0.28763186 

0.32303160 

0.35893158 


0.39535650 

0.43230305 

0.46973339 

0.50756428 

0.54565387 

0.58379085 

0.62169896 

0.65908255 





After line 47, the region of figure 11 is normal, and values at IT = -97 are desired. 
Since IT1 = -107 for IM = 48 to 50, the 11th value in these rows corresponds to IT = -97. 

Boundary 2 is easier for this example than boundary 1. At all values of IM, u is 
desired for IT = -90. At IM = 44, IT = -90 corresponds to the first value in the row, 
which is 0.01517. For IM = 50, we need the 18th value (-107 + 17 = -90), which is 

0. 21819. 

On boundary 3, values of u are desired from mesh line 44. For the first section 
of the boundary, values must once again be obtained from the periodic blade above. 

IT = -97 to -93 is the desired range. Adding NBBI = 28 gives IT = -69 to -65. The values 
corresponding to these IT's are circled in table VII (0. 79106 to 0. 97554). Subtracting 

1. 0 from each of these gives -0. 20894 to -0. 02446. For the final point on the blade sur- 
face (point a), u = 0 is used. 

The upper part of boundary 3 has two points. The first is again u = 0. The second, 
for IM = 44, is IT = -90. This point is the first stream-function value given for line 44, 
which is 0. 01517. 

Boundary 4 corresponds to IM = 50. Values are required from IT = -97 to IT = -90. 
Since IT1 = -107, the 11th to 18th values are desired. These are circled in table VII 
(-0. 02076 to 0.21819). 

After the final set of boundary values is given for boundary 4, a final data card must 
be given with zero for KBDRY and NSP. This signals the end of the data for MAGNFY. 


Output 

Generally, the MAGNFY output is similar to the 2DCP, TURBLE, or TANDEM out- 
put, but for the finer mesh. In MAGNFY the vertical mesh lines are numbered with 
IM = 1 for the left boundary to IM = MMM for the right boundary. The horizontal mesh 
lines are numbered with IT = 1 for the lower boundary to IT = ITMAX for the upper 
boundary. 

Sample output is given in table VTU for the first example. Since the complete output 
would be lengthy and is similar to that for 2DCP, TURBLE, or TANDEM, the only output 
reproduced here is that which differs from the output of these three programs. The main 
part of table VTU is the stream-function values (UBV) and values of p times the com- 
ponent of W normal to the boundary (RWBV) along the vertical and horizontal bound- 
aries for the finer mesh. This main part is followed by a table of calculated program 
constants. The variable names are all defined in the section Main Dictionary. 
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TABLE VHI. - EXAMPLE OF MAGNFY OUTPUT 


STREAM FUNCTION AND RHO*W-SUB-THETA ON HORIZONTAL BOUNDARIES 


LOWER HORIZONTAL BOUNDARY 


UPPER HORIZONTAL BOUNDARY 


M 

UBV 

RW8V 

UBV 

RWBV 

0.96607E— 01 

-0.20894 

-25.8414 

0.151 70E— 01 

-25.3671 

0 . 97203E-01 

-0.20220 

-25.9841 

0.218 10E—0 1 

-25.6806 

0.977 99E-01 

-0. 19550 

-26.0991 

0.28460E-01 

-26.0027 

0.98395E-01 

-0. 18884 

-26.1856 

0.351 19E— 01 

-26.3336 

0.98992E-01 

-0.18225 

-26.2423 

0.41788E— 01 

-26.6735 

STREAM FUNCTION AND RH0»W- 

-SUB-M ON 

VERTICAL BOUNDARIES 



LEFT VERTICAL 

BOUNDARY 

RIGHT VERTICAL 

BOUNDARY 

THETA 

UBV 

RwBV 

UBV 

RWBV 

-2.72084 

-0.20894 

34.8187 

— 0 . 20/70E— 01 

30.8383 

-2.71523 

-0.20090 

35.1925 

—0 • 1 3689E-0 1 

30.8007 

-2.70962 

-0.19276 

35.6345 

— 0.66120E— 02 

30.7562 

-2.70401 

-0.18451 

36.1446 

0.45354E-03 

30.7049 

-2.69840 

-0.17613 

36. 7226 

0. 75065E— 02 

30.6467 


CALCULATED PROGRAM CONSTANTS 


PITCH 

0.7853982 

MB 1 1 
-79 


HT 

0.5609987E-02 

MBOO 
16 


HM1 

0.5434444E-03 

MMM I TM AX 

31 36 


HM2 

0.5962105E-03 


HM3 

0.5962105E— 03 


NUMBER OF INTERIOR MESH POINTS = 870 


SURFACE 

SURFACE 

1 

2 

3 

4 


BOUNDARY VALUES 
BV 

0 . 

0 . 

-0.44431 

-0.44431 


Error Conditions 

The error conditions are as follows: 

(1) SPLINT USED FOR EXTRAPOLATION 
EXTRAPOLATED VALUE -X.XXX 

SPLINT is normally used for interpolation, but may be used for extrapolation in some 
cases. When this occurs, the above message is printed, as well as the input and output 
of SPLINT. Calculations proceed normally after this printout. 

(2) BLCD CALL NO. XX 

M COORDINATE IS NOT WITHIN BLADE 

This message is printed by subroutine BLCD if the m- coordinate given this subroutine 
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as input is not within the bounds of the blade surface for which BLCD is called. The value 
of m and the blade-surface number are also printed when this happens. This condition 
may be caused by an error in the integer input items for the program. 

The location of the error in the main program is given by means of BLCD CALL 
NO. XX, which corresponds to locations noted by comment cards at each MHORIZ, 

ROOT , and BLCD call in the program . 

(3) ROOT CALL NO. XX 

ROOT HAS FAILED TO CONVERGE IN 1000 ITERATIONS 
This message is printed by subroutine ROOT if a root cannot be located. The input to 
ROOT is also printed. The user should thoroughly check the input to the main program. 

The location of the error in the main program is given by means of ROOT CALL 
NO. XX, which corresponds to locations noted by comment cards at each MHORIZ and 
ROOT call in the program. 

(4) DENSTY CALL NO. XX 

NER(l) = XX 

RHO*W IS X.XXXX TIMES THE MAXIMUM VALUE FOR RHO*W 
This message is printed if the value of pW at some mesh point is so large that there is 
no solution for the values of p and W. This indicates a locally supersonic condition, 
which can be eliminated by decreasing WTFL in the original 2DCP, TURBLE, or 
TANDEM run to obtain new input boundary values for MAGNFY. 

If RHO*W is too large, MAGNFY still attempts to calculate a solution. This often 
permits an approximate solution to be obtained which is valid at all the subsonic points 
in the region. In other cases, the value of W is reduced at some of the points in ques- 
tion during later iterations, resulting in a valid final solution for these points. The pro- 
gram counts the number of times supersonic flow has been located at any point during 
a given run (NER(l)). When NER(l) = 50, the program is stopped. 

The location of the error in the main program is given by means of DENSTY CALL 
NO. XX, which corresponds to locations noted by comments cards at each DENSTY call 
in the program. 

(5) THE USER HAS FAILED TO SPECIFY WHICH TYPE OF INPUT HE IS USING 
The first card of input after the title card specifies whether a single or tandem blade is 
being considered. There must be a 1 or a 2 in column 5 of this card. 

(6) MMM GT 100, OR ITMAX GT 50 

This is printed if MMM > 100 or if ITMAX > 50. In this case either MAGFAC should be 
reduced, or a smaller region chosen. 

(7) ONE OF THE MH ARRAYS IS TOO LARGE 

This is printed if there are more than 100 intersections of horizontal mesh lines with any 
blade surface. In this case MAGFAC should be reduced, or a smaller region chosen. 
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(8) THE NUMBER OF INTERIOR MESH POINTS EXCEEDS 2000 

This is printed if there are more than the allowable number of finite-difference grid 
points. Either MAGFAC must be reduced, or a smaller region must be chosen. 

(9) SEARCH CANNOT FIND M IN THE MH ARRAY 

H this is printed, the value of m and the blade-surface number are also printed. The 
user should thoroughly check the input to the main program. 


PROGRAM PROCEDURE 

The program is segmented into seven main parts - the subroutines INPUT, PRECAL, 
COEF, SOR, SLAX, TANG, and VELOCY, called by the main program MAGNFY. In 
addition there are several other subroutines. All the subroutines and their relation are 
shown in figure 12. All information which must be transmitted between the seven main 
subroutines is placed in COMMON. 


I II.V.H | I INPUT 


IBDVINTh 

I splint! 


1PRECAU 




MAGNFY 


COEF 


SOR 


IMHORIZI ICOEFBBj - 


p fHRBT lAAKl 
IBDRY12 1 IBDRY341 


SLAX 


IWRITUl ISLAVBBI 


LEO 


ROOT 


I 


BLCD 


IS PLINE I 


IDENSTYF 


TANG] 


IVELOCYl 


ISEARCHl 


IvelbbI 


IPLOTMYI 

IPISTUGI 


|SPI^22| 

Figure 12. - Calling relation of subroutines. 


Although most subroutines have been changed from those in TANDEM (ref. 3), the 
subroutine names have not been changed. Even with the rather extensive revisions made 
in some subroutines, the general descriptions of the subroutines in reference 3 still apply 
with minor differences, except for INPUT and PRECAL. Therefore, the only subroutines 
described here are INPUT and PRECAL, plus two new subroutines, WRITU and BDVTNT. 

Most of the subroutines in MAGNFY use the same set of variables. And most of 
these variables are the same as those used in TANDEM. These variables are all defined 
in the Main Dictionary (p. 32). The subroutines using these variables are described prior 
to the Main Dictionary or in reference 3. The remaining subroutines are described after 
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MHORIZ SEARCH 

2555,8, 2326 (8) 


Figure 13. - Arrangement for overlay, showing octal storage requirements. 


the Main Dictionary or in reference 3, and variables are defined with each subroutine. 

The program can handle as many as 2000 mesh points on the IBM 2-7094-7044 direct- 
coupled system with a 32 768-word core. To provide for the handling of 2000 mesh points 
an overlay arrangement is used, as shown in figure 13. All subroutines not shown are in 
the main link. The total program storage requirement is 74044^ of which 53364^ is in 
COMMON blocks which are stored in the main link. The system storage requirement for 
our computer is 2764^ and unused storage is 750^. If there is a storage problem on 
the user's computer, the maximum number of mesh points should be reduced. The fol- 
lowing program changes are required to change the maximum number of mesh points: 

(1) Change the dimension of A, U, K, and RHO in the COMMON/AUKRHO statement 
to the maximum allowable number of mesh points. This statement occurs in most sub- 
routines . 

(2) In subroutine INPUT change the number of values of K and RHO to be initialized 
(the bound on the DO loop near statement 240). 

(3) In subroutine PRECAL change the statement following statement 210 and format 
statement 1130 to reflect the maximum allowable number of mesh points. The statement 
following statement 210 will cause the program to stop if there are too many mesh points. 

(4) Change the dimensions of W, RWM, and BETA in SLAX, SLAVBB, TANG, 
VELOCY, and VELBB. 

(5) If the number of mesh points is reduced below 1600, the EQUIVALENCE state- 
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ments in SLAX, SLAVBB, TANG, VELOCY, and VELBB must be changed. 

The first segment of the program is INPUT. This subroutine reads all input data 
cards, calculates constants, and initializes arrays. It also uses SPLINE interpolation 
on the input boundary values to obtain boundary values on the fine mesh. The next sub- 
routine is PRECAL, which calculates all quantities that remain constant for a single 
problem. INPUT and PRECAL are each called once for a given problem. The remaining 
subroutines are each called once for each outer iteration. The subroutine COEF calcu- 
lates the entries of the matrix A and the vector k of equation (A7) (ref. 1). These coef- 
ficients must be recalculated for each outer iteration. On the first outer iteration sub- 
routine SOR estimates an optimum overrelaxation parameter S2 on the first call if it is 
not given as input. The same value of £2 is used for each outer iteration. SOR then finds 
the linear solution to equation (A7) (ref. 1) with fixed coefficients by successive over- 
relaxation. Then subroutine SLAX calculates pW m - Subroutine TANG calculates pW^, 
and then pW and /3 throughout the region. Finally, the subroutine VELOCY calculates 
density p and velocity W throughout the region and on the blade surfaces, and plots the 
surface velocities. 


Conventions Used in Program 

In general, the same conventions are used in MAGNFY as were used in the TANDEM 
program (ref. 3). In addition the lower, upper, left, and right boundaries of the mag- 
nified region are numbered 1, 2, 3, and 4, respectively. Also, the lower and upper 
boundaries of the region must sometimes be considered as blade surfaces by the program. 
In these cases they are numbered 5 and 6, respectively. 


Labeled COMMON Blocks 

The labeled COMMON blocks are organized the same as for TANDEM (ref. 3), 
except for the omission of /SLA/. 


Subroutine INPUT 

Read and print first part of input . - The program first reads the input cards which 
are the same as those for TURBLE or TANDEM, or were modified from 2DCP. A 
description of this input is given in reference 1, 2, or 3. All the input data are printed 
as they are read in. 

Fill in dummy second blade for single-blade case. - When there is only one blade, 
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the arrays for the second blade (surfaces 3 and 4) are filled in with the data for the first 
blade. 

Ca lculate l arge mesh arrays . - The large mesh spacing and MV array are calculated 
first so that the RMI and RMO arrays can be calculated. 

T rans fer blade coordinates for tandem blade when necessary. - If a leading-edge 
region is being analyzed, it is always considered to be on the rear blade. Hence, if it 
is the leading edge of the first blade of a tandem or slotted blade, the arrays for the rear 
blade are filled with the data for the front blade. Similar considerations hold for a 
trailing-edge region, since it is always considered to be on the front blade. 

Calc ulate constants . - Geometrical and miscellaneous constants are calculated. 

Cal culate fin e-m esh m- coordinates . - The final MV array for the fine mesh can now 
be calculated. 

Read boundary values and interpolate for fine mesh . - Each set of coarse-mesh 
boundary value data is read. Then interpolated stream-function values for the fine mesh 
are calculated and printed. 


Subroutine PRECAL 

The calculation of X and other constants in PRECAL is no longer necessary in 
MAGNFY, since X is given as input and the other constants are not used in the calcula- 
tion. The remainder of the description of PRECAL in reference 3 is still valid. 


Subroutine WRITU 

WRITU prints the value of the stream function along a given vertical line between 
blades or boundaries. A label is printed with the value of IM for the mesh line and the 
value of IT(ITl) for the first printed stream-function value. 


Main Dictionary 

The Main Dictionary applies to the previously discussed subroutines. 

A array of coefficients of u (i.e. , elements of a^ of matrix A in 

eq. (A7) of ref. 1) 
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A12, A34 

AA 

AAA 

AANDK 

AATEMP 

ANS 

AR 

AZ 

B 

B12, B34 
BB 
BE 
BEH 

BESP 

BETA 

BETAH(BETAV) 

BETAI(BETAO) 

BETI(BETO) 


BLDAT 

BV 

BVTN 

BZ 


a 12’ a 34 i 11 e q # (A2) of ref. 1 
temporary variable in BLCD 
array used for temporary storage 
input variable 

temporary location for AANDK in SOR 

result of calls on ROOT in TANG and DENSTY in SLAVBB and 
VELBB 

input variable (from 2DCP, TURBLE, or TANDEM) 
aQ in eq. (A2) of ref. 1 

array containing stream- channel thickness b at the four points 
adjacent to a point for which AAK is called 

^12’ ^34 * n ec b (A2) re ^- 1 
temporary variable in BLCD 
array of values of b at vertical mesh lines 

array of values of b where horizontal mesh lines meet the four 
blade surfaces 

input variable (from 2DCP, TURBLE, or TANDEM) 
array of values of /3 at interior mesh points 

array of values of /3 where horizontal (vertical) mesh lines meet 
the four blade surfaces 

input variable (from 2DCP, TURBLE, or TANDEM) 

array of angles at tangent points of leading- (trailing-) edge radii 
with the four blade surfaces (see input BETH, 2, 3, 4 and 
BETOl, 2, 3, 4) 

input variable 

array of stream -function boundary values on the four blade 
surfaces 

input variable 

stream- channel thickness bp at point for which AAK is called 
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CDMBIT(CDMBOT) 

CHANGE 

CHORD 

CMM 

\ 

CP 

CPTIP 

DBDM 

DIST 

DMLR 


DTDM 

DTDMH(DTDMV) 

DTLR 

DUDM 

DUDT 

EM 

EMK, EMKM1 
ERROR 

ERSOR 

EXPON 

FIRST 

GAM 


temporary grid locations along meridional axis in INPUT 

change in value of stream function at a particular point during an 
iteration of SOR 

array containing the meridional chord distances of each of the 
four blade surfaces (see input CHORDF and CHORDR) 

temporary variable in BLCD 

specific heat at constant pressure, c 

2c T! 
p in 

array of slopes at vertical mesh lines of spline curve for stream- 
channel thickness 

meridional distance in SEARCH from a blade leading edge to 
where a horizontal mesh line meets a blade surface 

tolerance for mesh points near a boundary in the m-direction 
(If a mesh point is closer than DMLR to a boundary, the point 
is considered to be on the boundary. ) 

d0/dm along a blade surface in BLCD 

array of d0/dm where horizontal (vertical) mesh lines meet the 
four blade surfaces 

tolerance in 6 -direction (see DMLR) 

array of derivatives of stream function 3u/3m along horizontal 
mesh lines in meridional direction 

array of derivatives of stream function 3u/30 along vertical 
mesh lines in 8 -direction 

array of second derivatives of spline curves for each blade sur- 
face, calculated by SPLN22 in BLCD 

temporary variables for EM in BLCD 

maximum absolute value of change in u at any point for an over- 
relaxation (SOR) iteration 

input variable 

i/(r - i) 

initial value of some index 

input variable (from 2DCP, TURBLE, or TANDEM) 
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H 

HM1 

HM2 

HM3 

HT 

I 

IEND 

IH 


fflS 

IM 

IM1(IMT) 

IM2 

IMS 

IMSL 

IMSS 

IMTMl 

INF 

INIT 

INTVL 

IP 


array containing mesh spacing h between the point for which 
AAK is called and the four points adjacent to it 

mesh spacing in m-direction from upstream boundary through 
front blade 

mesh spacing in m-direction for overlapping portion of front and 
rear blades, or between blades for nonoverlapping case 

mesh spacing in m-direction through rear blade to downstream 
boundary 

mesh spacing in 0 -direction from blade to blade 

temporary integer variable in INPUT, PRECAL, SLAX, and 
SEARCH 

integer variable set equal to 1 when final convergence to a solu- 
tion is reached in the outer iterations on a given set of data 

array containing current number of intersections of horizontal 
mesh lines with each of the four blade surfaces as intersec- 
tions are located 

integer variable in BDRY34 and TANG for counting intersec- 
tions of horizontal mesh lines with blade surfaces 

index of mesh line in the meridional direction (m-direction) 

integer variable in TANG indicating vertical mesh line index of 
the first (final) point in region of a horizontal mesh line 

IMl + 1 

array containing total number of intersections of horizontal mesh 
lines with each of the four blade surfaces 

temporary variable in PRECAL 

temporary variable in PRECAL, VELOCY, and VELBB 

IMT - 1 

variable in PRECAL indicating when an infinite slope is located 
at a blade leading or trailing edge in a call on BLCD 

array used to indicate whether BLCD has been called previously 
on a given blade surface 

input variable 

index of mesh point 
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IP1,IP2,IP3,IP4 


IPL(IPU) 

IPLM1(IPUP1) 

IS 

IT 

ITF(ITL) 

ITER 

ITMAX 

ITOR 

ITV 


ITVL(ITVU) 

ITVMl(ITVPl) 

IV 

J 

K 

KA 


value of IP at the four points adjacent to the mesh point under 
consideration 

value of IP where a vertical mesh line meets a lower (upper) 
surface or boundary 

value of IP on a vertical mesh line adjacent to a lower (upper) 
surface in VELBB 

integer variable in SEARCH for indicating where a horizontal 
mesh line intersects a blade surface 

index of mesh line in 9 - direction 

input variable 

outer iteration counter 

maximum value of IT in magnified mesh region 

value of IT at origin of coordinates at leading edge of front blade 

array of horizontal mesh line indexes (IT) corresponding to 
itersections of vertical mesh lines with blade surfaces 
(ITV(IM, SURF) is the IT value for the mesh point in the region 
on vertical mesh line IM which is closest to blade surface 
SURF. If ITV s 1, the value is adjusted to 1 for a lower sur- 
face or 2 for an upper surface. If ITV ^ ITMAX, the value is 
adjusted to ITMAX - 1 for a lower surface or ITMAX for an 
upper surface. If a vertical line does not intersect a blade sur- 
face, its value of ITV(IM, SURF) is equal to -10 000. For the 
lower boundary (SURF = 5), ITV = 2; and for the upper boundary 
(SURF = 6), ITV = ITMAX - 1.) 

ITV of lower (upper) blade surface on a given vertical mesh line 

ITV of a blade surface in COEFBB for vertical mesh line to left 
(right) of the line under consideration 

array containing the value of IP at the base of each vertical mesh 
line 

temporary integer variable in INPUT 

array of constants; the vector k in eq. (A7) of ref. 1 

integer array indicating which of the four points surrounding a 
mesh point lie on a boundary 
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KAK 

KBDRY 

KK 

KKK 

LAMBDA 

LAST 

LER 

LMAX 

LOWER 

M 

MAGFAC 

MB1,MB2 

MBDYF, MBDYL 

MBI 

MBI2 

MBII 


MBIIM1 

MBIIP1 

MBIT, MBOT 

MBO 

MB02 

MBOO 


MBOOM1 

MBOOP1 

MH 


real array giving boundary values of points surrounding a mesh 
point next to a boundary 

input variable 

integer counter in BLCD 

array containing information used in plotting subroutine PLOTMY 

input variable 

final value of some index 

array indicating location of error messages printed by program 
maximum value of u^ n+ ^/u^ n for eq. (B2) of ref. 7 
integer variable representing one of lower blade surfaces, 2 or 4 
meridional coordinate, meters 
input variable 

temporary vertical grid line locations along meridional axis 
input variable 

input variable (from 2DCP, TURBLE, or TANDEM) 
input variable (from TANDEM) 

number of vertical mesh lines from left boundary to leading edge 
of a blade in the magnified region (If region surrounds trailing 
edge of rear blade of a tandem blade, MBII = 1000. ) 

MBII - 1 

MBII + 1 

temporary grid locations along meriodional axis 
input variable (from 2DCP, TURBLE, or TANDEM) 
input variable (from TANDEM) 

number of vertical lines from left boundary to trailing edge of a 
blade in the magnified region (If region surrounds leading edge 
of front blade of a tandem blade, MBOO = -1000. ) 

MBOO - 1 

MBOO + 1 

array of m- coordinates of intersections of horizontal mesh lines 
with the four blade surfaces 
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MLE 

MM 

MMM 

MMMM1 

MMLE 

MMMSP 

MR 

MRTS 

MSP 

MSPMM 

MV 

NBBI 

NBL 

NER 

NIP 
NOBL 
NP1, NP2 

NRSP 
NSP 
NS PI 

NSPM1 

OMEGA 

ORF 

ORFOPT 


array of m- coordinates of leading edges of the four blade sur- 
faces (see input MLE2) 

input variable (from 2DCP, TURBLE, or TANDEM) 

number of vertical mesh lines from the left to right boundaries 
of the magnified region 

MMM - 1 

temporary meridional distance in BLCD 

temporary meridional distance in BLCD 

input variable (from 2DCP, TURBLE, or TANDEM) 

integer switch in PRECAL indicating when infinite derivatives 
would be encountered in a call on MHORIZ 

array of m- coordinates of spline points for each blade surface 
measured from its leading edge (see input MSP1, 2, 3, 4) 

temporary meridional distance in BLCD 

array of m- coordinates of vertical mesh lines 

input variable (from 2DCP, TURBLE, or TANDEM) 

number of blades 

array indicating number of times certain error messages are 
printed by program 

number of interior mesh points 

input variable 

integer counters in VELOCY indicating number of plotted blade- 
surface velocities 

input variable (from 2DCP, TURBLE, or TANDEM) 
input variable 

array containing number of spline points on each of the four blade 
surfaces (see input SPLNOl, 2, 3, 4) 

NSP - 1 

input variable (from 2DCP, TURBLE, or TANDEM) 
input variable (from 2DCP, TURBLE, or TANDEM) 
upper bound for optimum £2 from eqs. (Bl) and (B2) of ref. 7 
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ORFTEM 

P 

PITCH 

R 

RATIO 

RBV 

RELER 

RHO 

RHOB 

RHOHB 

RHOIP 

RHOVB 

RI(RO) 

RM 

RMDTL2(RMDTU2) 

RMH 

RMI(RMO) 

RMM 

RMSP 

RW 

RWBV 


temporary storage for ORFOPT 

array containing information used in the plotting subroutine, 
PLOTMY 

6 - coordinate from blade to blade, 27r/NBL 

array of densities p at the four points adjacent to a point for 
which AAK is called 

value of u“ +1 /u“ for use in eqs. (B2) and (B3) of ref. 7 

array of densities p on the four boundaries of the magnified 
region 

maximum relative change in density at surface mesh points 
between two outer iterations 

array of densities p at interior mesh points 

temporary storage in VELBB for a value of p on a blade surface 

array of densities p at horizontal mesh line intersections with 
the four blade surfaces 

input variable (from 2DCP, TURBLE, or TANDEM) 

array of densities p at vertical mesh line intersections with the 
four blade surfaces 

array of leading- (trailing-) edge radii on the four blade surfaces 
(see input RI1, 2, 3, 4 and ROl, 2, 3, 4) 

array of r- coordinates of the mean stream surface radii at 
vertical mesh lines 

Jr(d0/dm)p at vertical mesh line intersections on lower (upper) 
blade surfaces 

array of r- coordinates of the stream surface radii where hori- 
zontal mesh lines meet the four blade surfaces 

array of r- coordinates of mean stream surface radii at the inlet 
(outlet) of the four blade surfaces 

temporary meridional distance in BLCD 

input variable (from 2DCP, TURBLE, or TANDEM) 

value of pW of a mesh point 

array of p times the velocity component normal to the boundary 
of the magnified region 
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RWM 

array of pW m where vertical mesh lines intersect the four 
blade surfaces 

RWT 

array of pW^ where horizontal mesh lines intersect the four 
blade surfaces 

RWMBV 

value of pW m at a mesh point along the lower or upper boundary 
of the magnified region 

RZ 

density Pq at point for which AAK is called 

S 

meridional distance between two adjacent blade-surface spline 
points in BLCD 

SI (ST) 

blade-surface number at the beginning (end) of a horizontal mesh 
line in TANG 

SAL 

array of values of sin a = dr/dm at each vertical mesh line 

SIGN 

integer constant in BLCD 

SPLNO 

number of input spline points on a blade surface 

SPM 

array of m- coordinates along a horizontal mesh line in TANG 

SRW 

integer code variable that will cause certain subroutines to write 
out useful data for debugging: 

SRW = 13, SPLINE will write input and output data. 

SRW = 16, SPLINT will write input and output data. 

SRW = 18, SPLN22 will write input and output data. 

SRW = 21, ROOT will write input and successive estimates of 

the root to which it is converging. 

STGR 

array of 0-coordinates of center of each trailing-edge radius with 
respect to center of its leading-edge radius (see input STGRF 
and STGRR) 

STRFN 

input variable 

SURF 

integer variable referring to each of the four blade surfaces, the 
lower boundary (5), or the upper boundary (6) 

SURVL 

input variable 

T1,T2 

elapsed time in clock pulses (1/60 sec) 

TBI 

tan 0 le 

TBO 

tan^te 

TGROG 

2yR/(y + 1) 
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TH 


THETA 
THK, THKMl 
THLE 

THSP 

TIME 

TIP 

TSP 

TV 

TWL 

TWLMR 

TWW 

U 

UBV 

UBVIN 

UNEW 

UPPER 

USP 

W 

WCR 


array in INPUT; also single variable in PRECAL and TANG con- 
taining 0-coordinates from leading edge of front blade to a 
horizontal mesh line 

0- coordinate of a point along a blade surface in BLCD 

temporary variables in BLCD 

array of 0- coordinates from origin of front blade to leading edge 
of each blade surface (see input THLE2) 

array of 0- coordinates of spline points for each blade surface 
measured from its leading edge (see input, THSP1, 2, 3, 4) 

elapsed time in minutes 

input variable (from 2DCP, TURBLE, or TANDEM) 

array of 0- coordinates of points along a vertical mesh line in 
SLAVBB 

array of 0- coordinates where vertical mesh lines meet the four 
blade surfaces 

2coA 

2u)A - (wr)^ 

2a>/w 

array of stream-function values at each mesh point, or of eigen- 
vector associated with calculation of ORFOPT 

array of values of stream function u at mesh points on boundary 
of magnified region 

input variable 

new value of stream-function estimate at a single point calculated 
by eq. (7) of ref. 1 

integer variable representing one of the upper blade surfaces, 

1 or 3 

array of values of stream function along a vertical or horizontal 
mesh line, including boundary points 

array of relative velocities W at unknown mesh points, also 
used for storing pW 

critical velocity on a blade surface 
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WMB 

WTB 

WTFL 

WTFLSP 

WWCRM 

WWCRT 

XDOWN 

YACROS 


array of pW m where vertical mesh lines intersect the four blade 
surfaces 

array of pWg where horizontal mesh lines intersect the four 
blade surfaces 

input variable (from 2DCP, TURBLE, or TANDEM) 
input variable (from TANDEM) 

array, ratio of blade-surface velocity (based on meridional com- 
ponents) to critical velocity 

array, ratio of blade-surface velocity (based on tangential com- 
ponents) to critical velocity 

array of m- coordinates where surface velocities are plotted 
array of surface velocities to be plotted 


Program Listing for Subroutines Using Main Dictionary 


COMMON SRW, ITER, I END , L FR ( 2 ) , NER ( 1) 

COMMON / AUKKHO/ A ( 20 0 0 , 4 ) , U ( 2000 ) , K { 2000 ) , RHO l 2000 ) 

COMMON /INP/ GAM.AR, TIP, PHOI P, WTFL, WTFLSP, 0MFGA, ORF, BETA I , BE TAG, 
1 N0BL » MB I , MRU, MB I 2 , MB02 , MM, NBB I , NBL , NR SP , MBDYF , MBD YL , I TF, l TL, 

2 BLOAT, A AN Dr , LRSOR.STRFN, INTVL, SURVL , MAGF AC , 

3MR( j 0 ) ,RMSP(f>0),BESP (BO) 

COMMON /CALCON/ MB 1 1 , M BOO, MM M , MB I I M 1 , M b I I P 1 , MaUOM 1 , MB00P 1 , MMMM 1 , 
1HM1 ,HM2,HM3,HT,DTLR, 0 MLR , P I TCH , CP , F XPON , TWW , CP T I P , TGRGG , TB I ,TBU, 
2LAMr DA, TWL , I TOR , I T MA X , N I P , I M S ( 4 ) , B V ( 4 ) , M V ( 1 00 ) , I V ( 1 0 1 ) , 
3UBV(100,4) ,I<WdV( 100,4) ,1 TV( 100,6) , I V ( 100,4) » DTDM V ( 100,4) , 

4BET-WI 100,4 ) ,MH( 100, 4 ) ,0T0MH( 100,4) , BE TAH( i 00,4) ,RMH( 100,4) , 

5BEH ( 100 ,4 ) ,RM ( 100 ) ,BF ( 10U) , DBDM ( 100 ) , SAL ( 100 ) , AAA { 100) 

COMMON /G FOMIN/ CHORD (4) , S TGR ( 4 ) , ML E ( 4 ) , THLE ( 4 ) , RM I ( 4 ) , RMU l 4 ) , 
1RI (4 ) , R0( 4 ) , BLT I ( 4) , RE TO (4 ) , NS PI ( 4 ) .MSP { SO, 4 ) , THSP ( 50,4) 

COMMON / R H 0 S / RHOHB { L 00 , 4 ) , R HO VB ( 1 00 , 4 ) , R B V ( 1 00 , 4 ) 

COMMON /BLCDCM/ EM ( 5 0 , 4 ) , I N I T ( 4 ) 

INTEGER BLDAT, AANDK, ERSOk, STRF N , SURVL , A A TEMP , SURF , 

1FIRST, UPPER, SI, ST, SRW 

REAL K,KAK,LAMI3DA,LMAX,MH,MLF,MR,MSP,MV,MVIM1 
CALL TIME I ( I l ) 

10 TEND = -1 
ITER = 0 
DO 20 SUR F= i ,4 
20 INI f ( SURF ) = 0 
CALL INPUT 
CALL PRECAL 
30 CALL COEF 
CALL SOR 
CALL TIMEKT2) 

TIME= ( T2-TD/3600. 

WRITEI6.1000) TIME 
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CALL SLAX 
CALL TANG 
CALL VELOCY 
CALL TIMEKT2) 

T I Mt= ( T2-T1 )/3600. 

WR[ TE(6, 1000) TIME 
IF ( I END ) 30,30,10 

1000 FORMAT (8HLTIME = ,F7 .4, 5H MIN.) 
END 


SUBROUTINE INPUT 

INPUT READS AND PRINTS ALL INPUT DATA CARDS AND CALCULATES HORIZONTAL 
MESH SPACING (MV ARRAY) 

COMMON SRW.ITER, IEND,LER(2) f NER( 1) 

COMMON / AUKKHG/ A ( 2000 , 4 ) , U ( 2000 ) , K ( 2000 ) , RHO { 2000 ) 

COMMON / I NP/ GAM, AR, T I P , RHO I P , WTFL , WTF L SP , OMEGA , DRF , BE TA I , BE T AO, 
1N0BL.M8I ,M8G,MBI2,MB02,MM,NBRI , NBL , NRSP , MBDYF , MBD YL , ITF, 1TL, 

2BLDAT , AANDK, ERSOR, STRFN, INTVL, SURVL, MAGFAC , 

3MR(50) .RMSP150) .BESP150) 

COMMON /CALC ON/ MSI I , M8Q0, MMM, MB 1 1 Ml , MB I IP 1 , M80UM1 ,MB00P1 ,MMMM1 , 
1HM1 ,HM2 , HM3 , HT , OTLR, DMLR ,P ITCH, CP, EXPON , TWW.CPTI P, TGROG, T8I , TBO, 

2 LAMBDA, TWL, I TOR , I TMAX , N I P, IMS(4),BV(4),MV(100),IV(101), 

3UB V ( 100,4 ) , RW6V ( 100, 4 ) , I TV ( 1 00, 6 ) , T V ( 100 , 4 ) , DTDM V ( 1 00 , 4 ) , 

48ETAVI 100 ,4) ,MH{ 100, 4 ) , D TDMH ( 100, 4 ) , BETAH{ 1 00, 4 ) , R MH { 100 , 4 ) , 

5BEH1 100,4) ,RM 1 100) ,BE ( 100) , DBDM ( 100) , SAL ( 100) , AAAI 100) 

COMMON /GEOMIN/ CHORD (4) ,STGR( 4) , MLE ( 4 ) , THLE ( 4 ) , RM 1(4), RMO ( 4 ) , 

1RI ( 4 ) , RO ( 4 ) ,BETI ( 4 ) , 8 ETO (4 ) , NSP I ( 4 ) , MSP ( 50 , 4 ) , TH SPt 50,4) 

COMMON /RHOS/ RHOHBt 1 00 , 4 ) , RHO VB ( 1 00 , 4 ) , RB V ( 100 , 4 ) 

DIMENSION BVIN( 100)»UBVIN( 100),TH( 100) 

INTEGER BLD AT, AANDK, ERSOR, STRFN , SURVL , A A TEMP , SURF , 
1FIRST,UPPER,S1,ST,SRW 

REAL K,KAK,LAMBDA,LMAX,MH,MLF,MR,MSP,MV,MVIM1 

RE AO A JD WRITE INPUT DATA 

WR I TE( 6 , 1000 ) 

READ (5,1100) 

WR I TE ( 6 , 1 100 ) 

READ (5,1010) NOBL 
WR I r E ( 6, 1110) NOBL 

IF (NOBL.EO. 1 .OR .NOBL . EQ .2 ) GO TO 10 
WRITE (6,1120) 

STOP 

C OLD TURBLE ( 2DCP ) OR TANDEM DATA 
10 IF (NOBL. EQ. 1) WRITE (6, 1130) 

IF (NOBL.EO. 2) WR I TE l 6 , 1 140 ) 

READ (5,1030) G AM , AR , T I P , RHO I P , WTFL , WTFL SP , OMEGA , URF 
WRITE (6, 1040) GAM, AR, TIP, RHO IP, WTFL, WTFL SP, OMEGA, ORF 
IF (NOBL.EQ.l) WRITE (6,1 150) 

IF (NOBL.EO. 2) WR I TE ( 6 , 1 160 ) 

READ (5 , 10 30) BET A I , BE T AO , CHORD ( 1 ) , S TGR ( 1 ) , CHORD ( 3 ) , S TGR ( 3 ) , 

1MLE ( 3 ) , THL E ( 3 ) 

WRI TE(6 , 1040 IBETAI , BE T AO .CHORD ( 1) ,STGR( 1 ) .CHORD ( 3) , STGRI 3) , 
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1MLE ( 3) ,THLE( 3) 

IF (N0BL.EQ.1) WRITE I 6, 1 170) 

IF (N0BL.E0.2) WR I TE { 6 , I 180 ) 

READ (5,1010) MBI ,MBO,MB 12,MB02,MM,NBBI,NBL,NRSP 
WR I TE ( 6 , L 0 1 0 ) MB I , MBO , M8 12 , M 802 , MM , NBB I , NBL , NRSP 
DO 20 J= 1 » 4 

IF (N03L. EQ. 1. AND. J. EQ.3 ) GO TO 30 
IF (J.EQ.l) WRITEI6, 1190) 

IF (J.E0.2) WRITE(6, 1200) 

IF (J.E0.3) WRITEI6, 1210) 

IF (J.EQ.4) WRITE(6, 1220) 

WRITE(6,1230) J,J,J,J,J 

READ (5,1030) R I ( J ) , RO ( J ) , BE TI ( J ) , BETO ( J ) , SPLNO 
WRI TE(6, 1040) RI ( J),RO( J >,8ETI ( J) ,8ET0( J ) , SPLNO 
NSP I ( J ) = SPLNO 
NSP = NSP I ( J) 

WRI TE(6, 1240) J 

READ (5,1030) ( MSP ( I , J ) , 1= l , NSP ) 

WRI IE (6, 1040) (MSP (I ,J), 1=1, NSP) 

WRI rE(6, 1250) J 

READ (5,1030) ( THS P ( I , J ) , I = l , NSP ) 

20 WR I TE ( 6 , 1 040 ) ( THSP ( I , J > , I = 1 , NSP ) 

30 WR I TE ( 6 , 1 260 ) 

READ (5,1030) (MR( I) , 1=1, NRSP) 

WRI TE<6, 1040) ( M R ( I) , 1 = 1, NRSP) 

WRI TE(6,1270) 

READ (5,1030) (RMSPl I ), 1=1, NRSP) 

WRI TE( 6, 104u ) (RMSP( I), 1=1, NRSP) 

WRI TE(6, 1280) 

READ (5,1030) (BESP( I ), 1=1, NRSP) 

WR I TE ( 6 , 1 040 ) (BESP( I ), 1=1, NRSP) 

C NEW MAGNFY DATA 
WRITE(6, 1290) 

READ (5,1010) 8LDAT, AANDK, ERSOR,STRFN, INTVL.SURVL 
WR I TE ( 6 , 1 020 ) BLOAT, A ANDK, ER SOR , STRF N , I NTVL , SUR VL 
WRITE(6,1300) 

READ (5,1010) MBDYF.MBDYL, ITF, ITL.MAGFAC 
WRITE (6, 1020) MBOYF, MBOYL, ITF, ITL.MAGFAC 
WR I TE( 6 , 1 310) 

READ (5,1030) LAMBDA 
WRI TE< 6, 1040 ) LAMBDA 
MLE(l) = 0. 

THLE(l) = 0. 

FOR SINGLE BLADE CASE, FILL IN OUMMY TANDEM BLADE 

IF (N0BL.EQ.2) GO TO 60 
WTFLSP= 0. 

CHORD! 3 ) = CHORD ( 1 ) 

STGR ( 3 ) = STGR(l) 

MLE ( 3 ) = 0. 

THL E ( 3 ) = 0. 

MB 12= MBI 
MB02= MBO 
RI ( 3)= RI ( 1 ) 

RI (4)= RI (2) 

RO ( 3 ) = RO ( 1 ) 

RO ( 4 ) = RO ( 2 ) 

BET I ( 3 ) = BETIIl) 
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BET I C ^ ) = BFTI(2) 

BET0(3)= BFTO(l) 

B E T □ ( 4 ) = BET012) 

NSP I (3 ) = NS P I ( 1 ) 

NSP I ( 4 ) = NS P I (2) 

NSP = NSP I ( 3) 

00 40 1=1, NSP 
MSP ( I » 3 ) = MSP! I ,1) 

40 THSP ( I » 3) = THSP ( I, 1) 

NSP= NSPI (4) 

DO 50 1=1, NSP 
MSP ( I » 4 ) = MSP ( I » 2 ) 

50 THSP ( I , 4 ) = THSP ( 1,2) 

CALCULATE LARGE MESH SPACING 

60 HM1 = CHORD! 1 ) /FLOAT ( M80-M8I ) 

IF (MB0-GT.MBI2 .AND. MBI.NE.MBI2) HM1= MLE ( 3 1 /FLOAT ( HBI 2-MB I ) 

HM2 = 1.E30 

IF (MBI2.NE.M60) HM2= ( CHORD (1) -ML E ( 3 )) /FLOAT ( MBO-MB I 2 ) 

HM3 = CHQRD(3)/FL0AT (MB02-MBI2) 

IF (MB0.GT.M8I2 .AND. MB0.NE.MBU2) HM3= ( CHORD { 3 ) +MLE l 3 ) -CHORD ( 1 ) ) 
1 /FLOAT ( MB02-MBC) 

P I TCH= 2.*3.1415927/FLOAr(NBL) 

HT = PITCH/FLOAT (NBBI ) 

CALCULATE LARGE MESH MV ARRAY, AND RM I , RMQ , AND BV ARRAYS 

MBOT = M l NO ( MBO , MB I 2) 

CDMBOT = AMI N1 (CHORD! 1 ) ,MLE(3 ) ) 

DO 70 I M= 1 , MBOT 
70 M V ( I M ) = FLOAT! IM-MBI )*HM1 
MB I T= MAX0(MBU,MBI2) 

CDMB I T= AMAaI (CHORD! 1 ) ,MLE( 3 ) ) 

DO BO I M=MBUT , MB IT 

BO MV ( I M ) = CDMbOT+FLOAT! IM-MBOT )*HM2 
DO 70 I M= MB L T , MM 

90 MV ( I M) = CDMBIT+FLOAT! IM-MBIT )*HM3 

CALL SPLINT ( MR, RMSP , NRS P, MV , MM, RM , S AL ) 

RM I ( 1 ) = RM ( MBI ) 

RM I ( 2 ) = R M ( M B I ) 

RM I ( 3 ) = R M ( M B I 2 ) 

RM I ( 4 ) = RM ( MB I 2 ) 

RMO ( 1 ) = RM(MBO) 

RMO ( 2 ) = RM(MBO) 

RMU ( 3 ) = RM(MB02) 

RMO ( 4) = RM ( MB02 ) 

BV! 1 )= 0. 

BV ( 2 ) = 0. 

B V ( 3 ) = -WTFLSP/WTFL 
BV ( 4 ) = BV ( 3 ) 

CALCULATE GEOMETRICAL CONSTANTS 

MB 11= ( MBI 2-MbDYF ) *MAGFAC+1 
MBOO= ( MBO-MBDYF )*MAGFAC+1 
M MM= (MBOYL-MdDYF) *MAGFAC+1 
I TMAX= ( I TL- I TF ) *M AGF AC+1 
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HM1 = HM1/FLUAT (MAGFAC ) 

HM2 = HM2/ FLOAT (MAGFAC ) 

HM3 = HM3/FL0AT ( MAGFAC > 

HT = HT/FLOAT (MAGFAC ) 

I TOK= -IT F+MAGFAC+ 1 
IF (NOBL.EQ.l) GO TO 130 

FOR TANDEM BLADE CASE, IF REGION SURROUNDS LEADING EDGE OF 

FRONT BLADE, STORF BLADE SURFACES 1 AND 2 INTO 3 AND 4 

IF (MBDYF.Gt.MBI .OR. MBDYL.LE.MBI ) GO TO 110 

IF ( ITF.GE.O.OR.ITL.LE.O) GO TO 110 

MB 11= (MB I-MBDYF ) *MAGFAC + 1 

MBOU= -1000 

HM2 = HM1 

HM3 = HM1 

CHURD ( 3 ) = CHORD ( 1 ) 

STGR ( 3 ) = STGR(l) 

MLE ( 3) = MLfc(l) 

THL E ( 3 ) = THLE(l) 

DO 100 J= 1 » 2 
R I ( J + 2 ) = R I ( J ) 

RO ( J + 2 ) = RU ( J ) 

BET I (J + 2) = BET I ( J) 

BETOU + 2) = BETO(J) 

RMI(J+2) = RMI(J) 

RMOU+2) = RMO(J) 

BV ( J+2 ) = BV ( J ) 

NSPKJ + 2) = NSPI(J) 

NS P = NS P I ( J ) 

DO 100 1=1, NSP 

MSP ( I , J+2 ) = MSP(I,J) 

100 THSP ( I , J + 2 ) = THSP ( I , J ) 

FOR TANDEM BLADE CASE, IF REGION SURROUNDS TRAILING EDGE OF 

REAR BLADE, STORE BLADE SURFACES 3 AND 4 INTO 1 AND 2 

110 IF (MBDYF.Gt.MB02.0R.MBDYL.LE.MB02) GO TO 130 

IF ( FLOAT ( I I F *MAGF AC ) *HT.GE. THL E ( 3 > + STGR (3).0R.FL0AT(I TL + MAGFAL) * 
1HT.LE.THLE(3)+STGR(3) ) GO TO 130 
MB 11= 1000 

MBOO= ( MB02-M6DYF )*MAGFAC+1 
HM1 = HM3 
HM2 = HM3 

CHORD ( 1 ) = CHORD ( 3 ) 

STGR ( 1 ) = STGR (3) 

MLE ( 1 > = MLE ( 3 ) 

THLE(l) = THL E ( 3 ) 

DO 120 J = 3 > 4 
RKJ-2) = R I ( J ) 

RO ( J-2 ) = RO ( J ) 

BET I { J — 2 ) = BETIIJ) 

BET 0 ( J — 2 ) = BE TO ( J ) 

RMI1J-2) = RMI(J) 

RMO ( J-2 ) = RMO(J) 

BVIJ-2) = BV ( J > 
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NSPHJ-2) = NSPI(J) 

NSP = MSP I ( J) 

DO 120 1=1, NSP 

MSP ( l, J-2 ) = MSP l I ,J ) 

120 THSPU.J-2) = THSP ( I , J ) 

C 

C CALCULATE MISCELLANEOUS CONSTANTS 

C 

130 CHORD! 2 ) = CHORD ( 1 ) 

CHORD! 4 ) = CHORD { 3 ) 

STGR ( 2 ) = STGR(l) 

STGR ( 4 ) = S TGR ( 3 ) 

MLE ( 2 ) = ML E ( 1 ) 

MLE ( 4 ) = ML E ( 3 ) 

THL E ( 2 ) = THLE(l) 

THL L ( 4 ) = THLEI3) 

DTLR= HT/ 1000. 

DMLR= AMI Nl (HM1 ,HM2, HM3) /1000. 

MB I I Ml = MBII-l 
MB I I PI = MB I 1 + 1 
MBOOM1 = MBOO-1 
M B 0 G P 1 = MBOtJ + 1 
M M M M 1 = MMM-1 
NEK (1) =0 
CP = AR/ (GAM-1. ) *GAM 
EXPON= 1. /(GAM-1.) 

TWW = 2. *OMEGA/WTFL 

CPT I P= 2 • *CP»T I P 

TGROG= 2.*GAM*AR/(GAM+1. ) 

C 

C CALCULATE FINE MESH MV ARRAY 

C 

M V ( 1 ) = MV(MbDYF) 

MBOT = MINOIMBI I.MBOO) 

MBO T = MAXO(MBUT, 1) 

DO 140 I M= 1 , MBOT 

140 M V ( I M ) = MV ( 1 ) +FLOAT( IM-1 )*HM 1 
MB I T = MAXO ( MB I I , MBOO ) 

MB I T= MINO(MBIT , MMM ) 

DO 150 IM=M60T,MBIT 

150 M V ( I M ) = MV ( MBOT ) +FLOAT ( I M-MBOT ) *HM2 
DO 160 I M=MB IT, MMM 

160 MV ( 1 M ) = MVtMBIT I+FLOAT ( I M- MBIT ) *HM 3 
DO 165 I M= 1 , MMM 

165 IF (ABS(MV(IM) J.LT.DMLR) MV(IM)=0. 

CALL SPLINT ( MR.RMSP* NRSP ,M V, MMM, RM , SAL ) 
CALL SPLINT(MR,BESP,NRSP,MV,MMM,3E,DBDM) 
C 

C FINISH READING NEW MAGNFY INPUT DATA 
C READ, COMPUTE, AND STORE BOUNDARY VALUES 
C 

DO 170 1=1,100 
DO 170 J= 1 , 4 
UBV(I.J) = 0. 

170 RWBV(I.J) = 0. 

180 READ (5,1010) K.BDRY ,NSP 
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IF (KBDRY. EU.O) GU TO 210 
WRITEI6, 1320) 

WR I TE ( 6 » 1020 ) KBDRY, NSP 
WRITE(6,1330) 

READ (5,1030) ( BVIN(I) , I =1 ,NSP ) 

WR I T E ( 6 , 1040 ) (BVIN( I ),I=1,NSP) 

WRI TE(6, 1340 ) 

READ (5,1030) ( UBV IN ( I ) , 1=1, NSP ) 

WRI TE ( 6 , 1 040 ) (UBVIN( I ), 1=1, NSP) 

IF (KBDRY. EQ. 3 .OR. KBDRY. EQ. 4) GO TO 190 
CALL BDVINT(BVIN,UBVIN»NSP,MV,KBDRY,DMLR»MMM) 

GO TO 180 

190 DO 200 I T = 1 , I TMAX 

200 TH { I T ) = FLOAT ( IT— ITOR )*HT 

CALL BDVINT(BVIN,UBVIN,NSP,TH,KBDRY,DTLR, I TMAX) 

GO TO 180 

210 DO 220 KBDRY = 1,2 
DO 220 I M= 1 , MMM 

220 RWBV(IM, KBDRY) = -RW B V ( I M, KB DR Y ) *WTFL/ B E ( I M > 

DO 230 KBDRY=3 , 4 
DO 230 I T= l , I T MAX 

230 RWBV (IT, KBDRY ) = RWBV ( I T .KBDRY ) *WTFL /BE ( I M ) /RM l I M ) 

IF(BLDAT.GT.O) WRITE (6,1350) (MV( IM ) ,UBV( IM, 1) ,RWBV( IM, 1 ) , 
1UBV(IM,2) ,RWBV(IM,2) , IM= l, MMM) 

I F ( BLOAT . GT . 0 ) WRITE (6,1360) ( TH( IT ) , UBV ( I T , 3 ) , RWBV ( I T , 3 ) , 

1 UBV ( I T , A ) ,RwBV( IT ,4) , I T= l » IT MAX) 

INITIALIZE ARRAYS 

DO 240 1=1,2000 
U ( I ) = 1. 

K ( I ) = 0. 

240 RHO ( I ) = RHOIP 

DO 250 I M = l , 100 
DO 250 SURF = 1,4 
RHOHBt IM,SURF)= RHOIP 
RHUVB( I M, SURF ) = RHOIP 
RBVUM.SURF) = RHOIP 
250 ITV(IM.SURF) = -10000 
DO 260 I M= l , 1U0 
I T V ( I M , 5 ) = 2 
260 ITVl I M , 6 ) = I TMAX— 1 

IF (MMM.LE.100.AND.ITMAX.LE. 100) RETURN 
WRITE (6,1370) 

STOP 

FORMAT STATEMENTS 


1000 FORMAT 
1010 FORMAT 
1020 FORMAT 
1030 FORMAT 
1040 FORMAT 
1100 FORMAT 
1 

1110 FORMAT 


( 1H1 ) 

(1615) 

(IX, 1617) 

(8F10.5) 

( 1X.8G16.7) 

( 80H 

(7X.4HN0BL/7X, 13) 


) 
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1120 FORMAT (29H1 NOBL HAS NOT BEEN SPECIFIED) 

1130 FORMAT ( 7 X , 3HGAM , 14X, 2hAR, 13X, 3HTIP, 12X.5HRH0IP, 12X,4HWTFL,2 /X, 
15H0MESA , 12X, 3H0RF ) 

1140 FORMAT (7X.3HGAM, 14X , 2HAR, 13X, 3HT I P , 1 2X , 5HRH0I P , 12X , 4HWTFL , 1 IX , 6HW 
ITFLSP, 10X,5n0MEGA, 12 X, 3H0RF) 

1150 FORMAT < faX , 3HBET A I , 10X , 5HBETA0, 1 1 X , 5HCHURD , 12X , 4HS TGR ) 

1160 FORMAT (6X,5HBETAI , 1 OX , 5HBET AO , 1 IX , 6HCH0R0F , 11X» 5HSTGRF, 10X, 
16HCH0RDR. 10X.5HSTGRR, 12X.4HMLER, l IX « 5HTHLER ) 

1170 FORM, AT ( 4X » 8HMB I MBO , 12X , 18HMM NBB I NBL NRSP) 

1180 FORMAT ( 4 1H MB I MBO MB12 MB02 MM NBB I NBL NRSP) 

1190 FORMAT (53HL BLADE SURFACE 1 — UPPER SURFACt - FRONT BLADE) 

1200 FORMAT ( 5 3HL BLADE SURFACE 2 — LOWER SURFACE - FRONT BLADE) 

1210 FORMAT (52HL BLADE SURFACE 3 — UPPER SURFACE - REAR BLADE) 

1220 FORMAT (52HL BLADE SURFACE 4 — LOWER SURFACE - REAR BLADE) 

1230 FORMAT I7X.2HRI, Il f 12X f 2HR0, II, 1 2X . 4HBE T I » I 1 , 1 IX , 4HBE TO, 1 1 , 1 1 X , 5HS 
1PLN0, I l ) 

1240 FORMAT { 7 X , 3HM SP , I 1 , 2X , 5HARR AY ) 

1250 FORMAT (7X,4HrHSP, II, 2X, 5H ARRAY ) 

1260 FORMAT ( 1 6HL MR ARRAY) 

1270 FORMAT ( 7X » 1 1HRMSP ARRAY) 

1280 FORMAT (7X,11H[>ESP ARRAY) 

1290 FORMAT (45HL 8LDAT AANOK ERSOR STRFN InTVL SURVL) 

1300 FORMAT (39HL MBDYF MBDYL ITF ITL MaGFAC) 

1310 FORMAT ( 7X,oH LAMBDA) 

1320 FORMAT (15HL KBDRY NSP) 

1330 FORMAT (7X.11HBVIN ARRAY) 

1340 FORMAT (7X.12HU8VIN ARRAY) 

1350 FORMAT ( 1 H 1 , 7X , 60HST R E AM FUNCTION AND RHO*W-SUb-THETA ON HORIZONTA 
1L BOUND AR IES//19X,55HLUWER HORIZONTAL BOUNDARY UPPER HORIZONTA 

2 L B'JUNDAR Y/dX, 1HM,6X,2(7X, 3HUBV, 12X, 4HRWBV.4X) / ( 1X.5G15.5) ) 

1360 FORMAT ( L H I , 7 X , 54HSTR E AM FUNCTION AND RHO*W-SUB-M ON VERTICAL BOUN 
1DAR IES//20X, 52HLEFT VCRTICAL BOUNDARY RIuHT VERTICAL BOUNDAR 

2Y/6X,5HTHETA,4X.2(7X, 3HUBV » 1 2X » 4HRWB V , 4X ) / ( 1 X, 5G 1 5 . 5 ) ) 

1370 FORMAT (28H1 MMM GT 100 CR ITMAX GT 100) 

FNU 


SUBROUTINE PRECAL 
C 

C PRECAL CALCUL A I Ed ALL REQUIRED FIXED CONSTANTS 
C 

COM, 'ION SR W, ITER, I END , L ER ( 2 ) , NER ( 1 ) 

COMMON / AUKRHU/ A ( 200 0 , 4 ) , U ( 2000 ) , K ( 2000 ) , RHO ( 2000 ) 

COMMON /IMP/ GAM, AR, T I P , RHO I P,WTFL,WTFLSP» OMEGA, OR F ,BETAI,BETAG, 
1 NOBL, MB I , MBO, M3 1 2, MBO 2, MM, NB BI , NBL, NRSP, MBDYF, MBDYL, ITF, ITL, 

2 BID AT, A AN DK, ERSOR, STRFN, INTVL , SURVL , MAGF AC , 

3 MR ( RO) , RMSP ( 50 ) , BESP { 50 ) 

COMMON /CAL CON/ MB I I , MBOG, MMM , M B I I M 1 , MB I I P 1 , MiiOQM 1 , MBOOP 1 , MM MM 1 , 
lHMl,HM2.HM3,Hr,DTLR, DMLR ,P I TCH , CP , EX PON , TWW , CP T I P , TGROG , TB I , TBO, 
2LAM-SDA, TWL, I TOR, ITMAX, NIP, IMS(4),DV(4),MV(L00),IV(L0l), 

BUB V ( 10 0,4 ) ,RWbV ( 100, 4 ) , I TV ( 1 00 , 6 ) , T V ( 1 00 , 4 ) , DTDM V ( 100,4) , 

4 BE T WI 100,4) ,MH( 100, 4 ) ,DTDMH( 100,4) , BET AH( 100,4) ,RMH{ 100,4) , 
5BEHI 100,4) ,RM( 100) ,BE ( 10 0) , D BDM ( 1 00 ) , S AL { 1 00 ) , AAA { 100) 

INTEGER B LO AT, AANLiK, ERSOR, STRFN, SURVL, AATEMP, SURF, 
lFIRdT,UPPFR,Sl ,ST, SRW 

REAL K,K AK, LAMBDA, LM AX, Ml!, ML E, MR, MSP, MV, MV I Ml 

EXTERNAL BL l , BL2 , BL3 , BL4 
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CALCULATE TV, ITV, IV, DTDMV AND BETAV ARRAYS 

TV, ITV, DTDMV AND BETAV ON BLADES 
SRW = 0 

MBIT = MAXOI I, (M8I-MBDYF )*MAGFAC+l ) 

MBOT = MINOIMMM.MBOO) 

IF (MBOT. LF. MBIT ) GO TO 20 
DO 10 I M=MB IT, MBOT 
LER ( 2 ) = 1 
C BLCO CALL NO. 1 

CALL BL1 (MVC IM) ,TV(IM, l) , DTDMV ( IM, 1 ) , INF ) 

ITV(IM,1)= I NT ( ( TV ( I M , 1 ) +DTLR ) /HT ) + I TOR 
IF( TV{ IM, 1 ) .GT.-DTLR )ITVIIM,1)=ITV(IM,1>+1 
ITV ( I M , 1 ) = MINOl ITV ( IM, 1) , ITMAX) 

I TV ( I M , 1 ) = MAXOI ITV I IM, l) ,2 ) 

BETA VI I M, 1 )=AT AN (DTDMV I IM, 1 )*RM( IM) >*57.295779 
LER I 2 ) = 2 

C BLCD CALL NO. 2 

CALL BL2IMVI IM)»TV(IM,2) .DTDMV! IM,2) , INF ) 

ITVI IM,2) =INT( (TV! I M , 2 ) -DTL R ) / HT ) + I TOR 
IF I TV! IM, 2 ) .LT .DTLR) I TV! IM,2 ) = 1 TV I I M , 2 ) - 1 
I TV I IM, 2 ) = MAXO I ITV I IM, 2) , 1 ) 

ITVI IM, 2 ) = MINOI ITV I IM, 2) , I TMAX-1 ) 

10 BETAVI IM, 2) = AT AN (DTDMV! IM,2)*RMI IM) )*57. 295779 
20 MBIT = MAXO I 1 , MB I I ) 

MBOT = MINOIMMM, (MB32-MBDYF ) *MAGFAC+1 ) 

IF I MB I T . GT . MBOT ) GO TO AO 
DO 30 I M=MB IT, MBOT 
LER I 2 ) = 3 

C BLCO CALL NO. 3 

CALL BL3IMVI IM) ,TV(IM,3) , DTDMV ( I M, 3 ) , I NF ) 

ITVI I M , 3 ) = INTI (TVIIM.3) +DTLR )/HT)+ITOR 
IF (TV! IM, 3 ) .GT.-DTLR ) I TV I I M , 3 ) =1 TV I I M, 3 ) +1 
I TV I I M , 3 ) = MINO I ITV I IM, 3) , ITMAX) 

I T V I I M , 3 ) = MAXOI ITVI IM, 3) ,2) 

BETAVI IM,3)= AT AN I DTDMV I IM,3)*RM(IM) )*57.295779 
LER 1 2) = A 

C BLCO CALL NO. 4 

CALL BL41 MV I IM) ,TV( IM,4> , DTDMV I IM,4> , INF) 

ITVI IM,4) = INTI I TV! I M, 4) -DTLR )/HT) + ITOR 
IF I TV! IM,4 ) .LT.DTLR) I TV! IM,4 ) = I TV l I M , 4 ) - 1 
ITVI I M , 4 ) = MAXOI ITVI IM, 4) , 1 ) 

ITVI I M , 4 ) = MINOI ITV ( IM, 4) , I TMAX- l ) 

30 BE T A V I I M , 4 ) = AT AN I DTDMV I IM , 4 ) *RM I I M > ) *57 . 295 77 9 
C IV ARRAY 

40 I V I 1 ) = 0 
I V I 2) = 1 

MBOT = MI NO { MB 1 1 M 1 , MBOO ) 

IFIMBOT.LT. 2) GO TO 60 
00 50 I M= 2, MBOT 

TV I IM+1) = I V I IM ) + IT V I IM,2 ) — I T V C IM, 1 )+ITMAX-l 
50 IF I ITV I IM, 1 ) .EQ.— 10000) I VI IM + i )= I V I IM + 1 )- 1 
60 MBIT = MAXO ( 2 , MB I I ) 

MBOT = MI NO (MMMM1, MBOO ) 

IF! MB IT. GT. MBOT )G0 TO 80 
DO 70 I M=MH IT, MBOT 

70 I V ( I M+ 1 ) = I V ( I M ) + ITV I IM, 4) + 1 TV I IM,2) — ITVI I M , 3 ) — I TV ( I M, 1 ) +ITMAX 
80 MBIT = MAX0(2,MB00P1 ) 

MBOT = MINO I MMMM1 , MB I IM1 ) 
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IF(MBIT.GT.MBUT) GO TO 100 
DO 90 I M= Mu IT, MBOT 
90 I V ( IM+1 )= IV( IMJ+ITMAX-2 
100 MBIT = MAX0(M6II,MB00P1) 

IFIMBIT.GT.MMMMl )G0 TO 120 
DO 110 I M=MB I T » MMMM1 

I V { IM+l ) = IV( I M ) + 1 T V( I M , 4 )- I T V ( I M, 3 ) + I TM AX- 1 
110 IF ( I T V C IM,3) .EQ.-10000) I V ( I M+ 1 } = I V { I M+ 1 ) - 1 
120 N I P= I V { MMM ) — 1 

WRITE (6,1020) PITCH, HT, HM 1 , HM2, HM3 
WRITE (6,1030) MBI I, MBOO.MMM, ITMAX.NIP 
WRITE (6,1040) (SURF, BV( SURF), SUR F= 1 , 4 ) 

IF (RLDAT.LF.O) GO TO 140 

MBIT = MAXO( 1, (MBI— MBDYF )*MAGFAC+1 ) 

MBOT = MI MO (MMM, MBOO ) 

WRITE ( 6 , 10 50) 

DO 130 SURF = 1,3,2 
I = SURF+1 

IF (MBIT.Lt.MbOT ) WR I T£ ( 6, 1060 ) SURF , I , { MV ( I M) , T V ( I M, SURF ) , DTDMV ( I 
1M, SURF ) , TV ( I M, I ) , DTDMV ( IM, I ) , I M=MB I T , MBOT ) 

MBIT = M A X 0 ( l.MBII ) 

130 MBOT = MI NO (MMM, ( MBQ2-MBDYF) *MAGFAC+ l ) 

WRI TE( 6, 10 70 ) ( IM,MV( I M) , RM ( IM) .SAL ( IM) ,Bh( IM) ,DBDM( IM) , IM = 1 , MMM) 
WRITE (6 , 1080) ( IM, IV ( IM) , ( ITV( I M, SURF) , SURF=1,4) , IM=1,MMM) 

C 

C CALCULATE MH AND DTDMH ARRAYS. 

C 

140 IMS ( 1 ) =0 
MRTS = 0 

MBIT = MAXO ( 1 , { MB I -MBDYF )*MAGFAC + 1 ) 

MBOT = MI NO ( MMM , MBOO ) 

LER ( 2 ) = 5 

C 8LC0 AMD ROOT (VIA MHORIZ) CALL NO. 6 

CALL MHORI Z (MV, I TV (1,1), RL 1 , MB I T , MBU T , I TOR , HT , D TLR , 0 , I MS ( 1 ) , 

1 MH ( 1,1 ) , DTDMH ( 1,1) »MR TS ) 

IF ( ITV ( MBOO, 1 >-ITV(MBOO, 2) .NF.2) GO TO 160 

I MSL = I M S ( 1 ) + l 

MH ( (MSL, 1 )= MV (MBOO) 

D T D i" H ( I MS L , 1 ) = -1.E10 
I MS ( 1 ) = I MSL 
150 I MS ( 2 ) =0 
MRTS = 0 
LER ( 2) = 6 

C BLCD AND ROUT (VIA MHORIZ) CALL NO. 6 

CALL MHOR IZ ( MV, ITV ( 1, 2) , bL2 , MB I T , MBOT, I TOR. HT, DTLR , 1 , I MS ( 2 ) , 

1MH ( 1,2) .DTDMH (1,2), MRTS) 

I MS ( 3 ) =0 

I F ( I TV ( MB I I,3)-ITV(MBI I, 4) . NE . 2 ) GO TO 160 
MRTS = 1 
I MS { 3 ) = 1 

MH( 1,3)= MV ! MB I I ) 

DTDMHl 1,3)=1. E10 
160 MBIT = MAXO ( 1 , MB I I ) 

MBOT = MI NO ( MMM, ( MB02-M8DYF ) *MAGFAC+1 ) 

LER ( 2 ) = 7 

C BLCD AND ROOT (VIA MHORIZ) CALL NO. 7 

CALL MHOR I Z ( MV , I TV (1,3), BL3 , MB I T , MBOT, I TOR , HT , DTLR , 0 , I MS ( 3 ) , 

1MH( 1,3) , DTDMH ( 1,3) ,MRTS) 

I MS ( 4 ) -0 

IF !ITV(MBII,3)— ITV(MBII,4).EQ.2> MRTS=1 
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LER ( 2 ) = 8 

C BLCD AND ROOT (VIA MHORIZ) CALL NO. 8 

CALL MHORIZ (MV, I TV ( 1 , 4 ) , BL4, MB I T , MBOT, I TOR , HT, DTLR , I, I MS (4) , 

IMHt 1,4) , DTDMH (1,4) ,MRTS ) 

I = MAXOl IMS!1),IMS(2),IMS(3),IMS(4) ) 

IF (I.LE.lOO) GO TO 170 
WR I TE ( 6 , 1090 ) I 
STOP 

CALCULATE RKH, BEH, AND BETAH ARRAYS 

170 IF (BLDAT.GT.O) WRITE (6, 1100 ) 

DO 190 SURF= 1,4 

CALL S PL I NT ( MR, RMS P, NRSP ,MH( 1, SURF ) , IMS { SURF ) , RMH ( l.SURF) , AAA ) 

CALL S PL I NT ( MR ,BESP»NRSP»MH( l.SURF ) , IMS ( SURF ) , BLH ( 1 , SURF ) , AAA) 

I MSS= I MS ( SURF ) 

IF (IMSS.LT.l) GO TO 190 
DO 180 IHS= 1 , I MSS 

180 BET AH ( I HS » SURF ) = AT AN ( DT DMH ( IHS , SURF ) *RMH( IHS.SURF) >*57.295779 
IF (BLDAT.GT.O) WR I T E ( 6, 11 10 ) SURF , ( MH ( I M , SURF ) , RMH ( I M, SURF ) , 

1BEH ( IM, SURF ) , BETAH ( IM, SURF ) , DTDMH ( IM, SURF ) , IM=1 , I MSS ) 

190 CONTINUE 

IF (BLDAT.LE.O) GO TO 210 
WRITE (6,1120) 

DO 200 I T= 1 , ITMAX 
TH = FLOAT ( IT-ITOR)*HT 
200 WRITE (6,1010) I T , TH 
210 WRI IE (6, 1000) 

IF (NIP. LE. 2000) RETURN 
WRITE (6,1130) 

STOP 

1000 FORMAT ( 1 H 1 ) 

1010 FORMAT (4X, I4.G16.5) 

1020 FORMAT ( l H 1 / / / / / 5X , 2 8 HCALCUL AT ED PROGRAM CON STANTS //5X , 5HP I TCH , 

1 13X.2HHT , 13X.3HHM1, 13X.3HHM2 , 1 3X , 3HHM 3/ 1 X , 5G16 .7 ) 

10 30 FORMAT ( / 5X , 4HMB I I , 10X,4HMB00, 10X, 3HMMM , 10X.5HITMAX/3X, 15, 9X, 15, 
19X, I 5,9X, I 5///5X, 33HNUMBER OF INTERIOR MESH PUINTS = ,15) 

1040 FORMAT (////// 5X , 23HS URF ACE BOUNDARY V AL UE S / / 5 X , 7H SURF AC E , 7 X , 2HB V/ 
1 (5X, I4,4X,F10.5) ) 

1050 FORMAT ( 1 H 1 , 6X , 62HBL A DE DATA AT INTERSECTIONS OF VERTICAL MESH LIN 
1ES WITH BLADES) 

1060 FORMAT ( 1 HL , 22X , 1 3HB L ADE SURFACE , I 2 , 15X , 1 3HBLADE SURF ACE , I 2/ 7X , 
11HM, 14X.2HTV, 11X.5HDTDMV »12X,2HTV» l IX , 5H0TDMV/ { 5G 15. 5 ) ) 

1070 FORMAT ( 1 H 1 , 1 3X , 44HS T RE AM SHEET COORDINATES ANu THICKNESS TABLL/ 
12X, 2HIM.7X, 1HM, 14X.1HR, 13X, 3HSAL, 13X, 1HB , 1 2X , 5HDB/DM/ 

2 ( IX , 13 , 5G 15.5) ) 

1080 FORMAT (4H1 IM.9X.8HIV A HR A Y , 32X , 9H I T V ARR AY/38X , 5HBLADE/3 7X , 7HSUR 
1FACE.3X, 1H1, 5X, 1H2,5X , 1H3, 5 X , 1H4/ 39X , 3HN0. / ( 1X,I3,5X,U0,25X, 

24 ( I 4 ,2X ) ) ) 

1090 FORMAT (35H1 ONE OF THE MH ARRAYS IS TOO LARGE / 7HL I T HAS, 15, 

18H POINTS) 

1100 FORMAT ( 67H1M COORDINATES OF INTERSECTIONS OF HORIZONTAL MESH LINE 
IS WITH BLADE) 

1110 FORMAT (25HLMH ARRAY - BLADE SURFACE , I2//15X.2HMH, 19X, 3HRMH , 19X, 

1 3HBEH, 1 8X » 5HBET AH , 17X , 5HDTDMH/ l 5G22. 4 ) ) 

1120 FORMAT (43H1THETA COORDINATES OF HORIZONTAL MESH L I NE S//6X , 2H I T , 

1 5X , 5HTHET A ) 

1130 FORMAT (48H THE NUMBER OF INTERIOR MESH POINTS EXCEEDS 2000) 

END 
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SUBROUTINE COEF 


COEF CALCULATES FINITE DIFFERENCE COEFFICIENTS, A, AND CONSTANTS, K, 
AT ALL UNKNOWN MESH POINTS FOR THE ENTIRE REGION 

COMMON SRW,ITER,IEND,LER(2),NER( 1) 

COMMON / AUKRHQ/ A ( 2000 , 4 ) , U( 2000 ) , K ( 2000 ) , RHOI 2000 ) 

COMMON /INP/ GAM,AR,TIP,RHOIP,WTFL, WTFLSP, OMEGA, ORF , BETAI , BETAO, 
1N0BL,MBI ,MB0,MBI2,MB02,MM,NBBI , NBL , NRSP ,MBDYF , MBDYL , I TF , ITL, 
2BLDAT, AANDK,ERSOR,STRFN, INTVL , SUR VL , MAGF AC , 

3MR( SO) ,RMSP(50) ,BESP( 50) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM , MB 1 1 Ml , MB I I P 1 , MB00M1 , MBQ0P1 ,MMMM1 , 
IHM1,HM2,HM3,HT,DTLR,DMLR ,P ITCH, CP, EXPON , TWW, CP TI P , TGROG, TB I , TBO, 
2LAMBDA,TWL,IT0R, ITMAX.NIP, IMS(4),BV(4) , MV (100) ,IVt 10 I) , 

3UBV( 100 ,4) ,RWBV( 100, 4 ) , I TV ( 100, 6 ) , T V ( 100, 4 ) , DTDMV ( LOO, 4 ) . 

4BETAVI 100,4) , MH ( 100, 4 ) , DTDMH ( 100, 4 ) , BET AH { 100,4) ,RMH( 100,4) , 

58EH( 100,4) ,RM( 100) , BE ( 100) , DBDM ( 100 ) , SAL ( 1 00 ) , AAA ( 100) 

COMMON /HRBAAK/ H ( 4) , R ( 4 ) , B ( 4 ) , KAK ( 4 ) , K A { 4 ) , I H { 4 ) , RZ ,BZ 
INTEGER BLOAT, AANDK, ERSOR, STRFN , SURVL , AATEMP , SURF , 

1 FIRST, UPPER, SI, ST, SR W 

REAL K, KAK, LAMBDA, LMAX,MH, ML E, MR, MSP, MV , MV I Ml 
C INITIALIZE ARRAYS 
ITER = ITER+1 

I H ( 1 ) = MAX0(0,ITV(2, 1)— ITV( 1,1) ) 

IH ( 2 ) = MAX0(0»ITV(1,2)-ITV{2,2) ) 

IH ( 3 ) = MAX0(0,ITV(2,3)— ITV( 1,3) ) 

I H ( 4 ) = M A XO ( 0 , I TV ( 1 , 4 ) — IT V ( 2, 4 ) ) 

I F ( ITV(MBII,3)-ITV(MBII,4).EQ.2) IH(3) = 1 
C INCOMPRESSIBLE CASE 

IF(GAM.NE. 1.5.0R.AR.NE.1000..0R.TIP.NE.1.E6) GO TO 10 
I END = 0 
GO ro 20 

C ADJUSTMENT OF PRINTING CONTROL VARIABLES 
10 IF ( ITER.NE. 1 .AND. ITER.NE.2 ) GO TO 20 
AANDK = AANDK— 1 
ERSOR = ERSOR— 1 
STRFN = STRFN— 1 
INTVL = INTVL-1 
SURVL = SURVL-1 
20 I F { IEND.NE.O) GO TO 30 
AANDK = AANDK+2 
ERSOR = ERSOR+2 
STRFN = STRFN+2 
INTVL = INTVL+2 
SURVL = SURVL+2 
C 

C CALL COEFBB THROUGHOUT THE REGION 
C 

C FRONT BLADE 

30 MBOT = MINOIMBI I Ml, MB 00) 

IF (MB0T.LT.2) GO TO 50 
DO 40 I M= 2 * MBUT 
CALL COEFBB (IM,5,2) 

40 CALL COEFBB ( IM, 1,6) 

C OVERLAP REGION 

50 MBIT = MAXO ( 2 , MB II) 

MBOT = MIN0(MMMM1,MB00) 

IF (MBIT. GT. MBOT) GO TO 70 
DO 60 I M=MB I T , MBOT 
CALL COEFBB (IM,5,4) 
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CALL COFFBB (IM,3,2) 

60 CALL COEFBB (IM,1,6) 

GO TO 90 

C NON-OVERLAP REGION 

70 MBIT = MAXO ( 2 » M BOO PI ) 

MBOT = MIN0(MBIIM1»MMMM1 ) 

IF (MBIT. GT. MBOT) GO TO 90 
DO BO I M= MB I T » MBOT 
80 CALL COEFBB (IM,5,6) 

C REAR BLADF 

90 MBIT = MAX0(MBII .MBOOPI) 

IF ( MB I T. GT . MMMM1 ) GO TO 110 
DO 100 I M=MB IT, MMMM1 
CALL COEFBB IIM,5,4) 

100 CALL COFFBB (IM,3,6) 

SPECIAL CASES - POINTS J OR C ARE MESH POINTS 

POINT J 

110 IF ( ITVIMBII »3)— ITV( MBI I ,4) . NE . 2 ) GO TO 120 
IT = I T V ( MB I I « 4 ) + 1 
IP = I PF ( MB I I M 1 , IT) 

K ( I P ) = K( IP)+A( IP, 4) *BV(4) 

A ( I P , 4 ) = 0. 

POINT C 

120 IF ( ITV(MBOO, 1 ) - I TV ( M BOO , 2) .NE.2) RETURN 
IT = I T V ( MBOO » 2 ) + l 
IP ■- I PF ( MBOOPI , IT) 

A ( IP. 3) = 0. 

RETURN 
END 


SUBROUTINE COEFBB ( IM , UPP ER , L OW ER ) 

COEFBB CALCULATES FINITE DIFFERENCE COEFFICIENTS, A, ANO CONSTANTS, K 
ALONG ALL VERTICAL MESH LINES WHICH INTERSECT BLADES 

COMMON / AUKRHG/ A ( 2000 , 4 ) , U ( 2000 ) , K ( 2000 ) , RHO ( 2000 ) 

COMMON /l UP/ GAM, AR, TIP, RHO I P , WT FL , WTFL SP , OMEGA , ORF , BE TA I , BE T AO, 
1N0BL.MBI , MBO.MBI 2,MB02,MM, NB B I , N8L , NR SP , MBD YF , MBD Y L , I TF , I TL, 

2 BLOAT, AANDK, ERSOR, STRFN, INTVL , SURVL , MAG FAC , 

3MR( 50) ,RMSP(50) , BESP ( 50) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM , MB I I M 1 , MB I I P 1 , MB00M1 , MBOOP 1 , MMMM 1 , 
1HM 1,HM2,HM3,HT,UTLR«DMLR,PITCH,CP» E XPON ,TWW»CPTIP* TGROG, TBI , TBO, 
2LAMBDA, TWL, ITOR, I TMA X , N I P, I M S ( 4 ) , BV ( 4 ) , MV l 1 00 ) , I V ( 10 1 ) , 

3UBV ( 100,4) , RWB V ( 100, 4 ) , I TV ( 1 00 , 6 ) , T V ( 1 00 , 4 ) ,DTOMV( 100,4) . 

4BETAV( 100,4) ,MH( 100,4) , D TDMH ( 1 00, 4), BET AH ( 100, 4), RMH (100, 4), 

5BEH( 100,4) ,KM( 100) ,BE ( 100) ,OBDM( 100) , SAL ( 100) , AAA ( 100) 

COMMON /HRHAAK/ H ( 4 ) , R ( 4 ) , B ( 4 ) , K AK ( 4 ) , K A ( 4 ) , I H ( 4 ) , R Z , B Z 
INTEGER BL DAT, AANDK, ERSOR, STRFN , SURVL , A A TEMP , SURF , 

IF I RST, UPPER, SI, ST, SR W 

REAL K,KAK,LAMBDA,LMAX,MH,MLF,MR,MSP,MV,MVIM1 
I F ( ITV( IM, UPPER). GT.ITVI IM, LOWER) ) RETURN 
I TVU= MAXO ( I TV( IM, UPPER) ,2 ) 

I T VL= M I NO ( I TV ( I M, LOWER I .ITMAX-l) 

IF ( ITVU.GT. ITVL) RETURN 
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I T= ITVU-1 

IPU = I PF ( IM, I TVU ) 

I PL = IPU+ITVL-ITVU 
00 80 I P= I PU » I PL 

it= ir+i 

CALL HRB ( IM, IT, IP) 

DO 10 1=1,4 
KAK ( I ) = 0. 

10 KA { I ) = 0 

FIX HRB VALUES FOR CASES WHERE MESH LINES INTERSECT BLADES 
OR BOUNDARIES 

IF (IT.NE.2) GO TO 20 
K AK ( 1 ) = UBV ( I M , 1 ) 

KA ( 1 ) = 1 

20 IF ( IT.NE. ITMAX-1) GO TO 30 
KAK ( 2 ) = UBV { IM, 2 ) 

KA ( 2 ) = 1 

30 IF ( IM.NE.2) GO TO 40 
KAK I 3 ) = UBV { IT, 3) 

K A ( 3 ) = 1 

40 IF ( IM.NE.MMMM1) GO TO 50 
KAK ( 4) = UBV { I T , 4 ) 

K A ( 4 ) = 1 

50 IF ( IT. EQ. ITVU- AND. UPPER .NE.5) CALL BDRY12 ( 1 , I M, I T , UPPER ) 
IF ( IT.EQ. ITVL.AN0.L0WER.NE.6) CALL BDR Y 12 t 2 , I M , I T , LOWER ) 
ITVM1= ITV( IM-1, UPPER ) 

I TVP 1= I T V { IM+l, UPPER ) 

IF { I T V ( IM, UPPER) .EQ. -10000) GO TO 55 
IF ( IT.LT. I IVMl ) CALL BDRY 34 ( 3 , I M, UPPER ) 

IF ( IT.LT. ITVP1) CALL BDRY 34 ( 4 , I M , UPPER ) 

55 l T V M 1 = ITVlIM-1, LOWER) 

I T VP 1= ITV( IM+l, LOWER ) 

IF ( ITVM1 .EU. -10000) GO TO 60 
IF (IM.EQ.MQII. AND. LOWER .EQ. 4 ) GO TO 60 
IF (IT.GT.ITVM1) CALL BDRY 34 ( 3 , I M , LOWER ) 

60 IF ( ITVP1.EQ. -10000) GO TO 70 

IF ( IM.EQ.MBOO.AND.LOWER.EQ.2) GO TO 70 
IF ( IT.GT. ITVP1) CALL BDRY34 ( 4 , I M , LOWER ) 

COMPUTE A AND K COEFFICIENTS 

70 CALL A AK ( IM,IP) 

DO 80 1=1,4 

K ( I P ) = K( IPJ+KAKU )*A( IP, I ) 

80 IF (KA(I).EQ.l) A(IP,I)=0. 

RETURN 

END 


54 


I 



noon onono 


SUBROUTINE HRBI IM, IT , IP) 

HRB CALCULATES MESH SPACING. H, DENSITIES, RZ AND R, AT GIVEN AND 
ADJACENT POINTS, AND STREAM SHEET THICKNESSES, BZ AND B, AT GIVEN 
AND ADJACENT POINTS 

COMMON / AUKRHG/ A ( 2000 , 4 ) , U ( 2000 ) , K ( 2000 ) , RHO l 2000 ) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM, MB 1 1 Ml , MB 1 1 P 1 , M6Q0M1 , MBOOP I ,MMMMI, 
1HM 1 ,HM2 ,HM3 , HT »DTLR, DMLR ,P ITCH, CP, EXPON , TWW, CP T I P , TGROG, TB I ,TBO» 
2LAMBDA,TWL, IT OR , IT MAX, NIP, IMS(4),BV(4),MV( 100) • I V( 101 ) , 

3UBV < 100,4 ) ,RWBV( 100,4), I TV ( 100, 6), TV { 100,4) ,DTDMV( 100,4) , 

4BET AVI 100,4) ,MH( 100,4 ) , DTDMH 1 1 00 , 4 ) , BET AH( 100, 4 ) , RMH ( 1 00 , 4 ) , 
5BEHI 100,4) ,RM( 100) , BE ( 100) ,DBDM( 100) , SAL ( 1 00 ) , AAA ( 100 ) 

COMMON /RHOS/ RHOHBI 1 00 , 4) , RHOVB ( 100 , 4 ) , RB V ( 100 , 4 ) 

COMMON /HRBAAK/ H ( 4) , R ( 4 ) , B { 4 ) ,KAK ( 4 ) , KA ( 4 ) , IH t 4 ) , RZ , BZ 
INTEGER BLOAT, AANOK, ERSOR, STRFN, SURVL , AATEMP , SURF , 

1 FIRST, UPPER, SI, ST, SRW 

REAL K,KAK, LAMBDA, LMAX.MH, ML E, MR , MSP .MV, MV I Ml 
H { 1 ) = HT*RM( IM) 


H ( 2 ) 

= 

HT*RM( IM) 




H ( 3 ) 

= 

M V ( I 

M) - MV ( 

IM-1 

) 


H 1 4 ) 

= 

MV( 

I M+ I ) —MV 

(IM) 



RZ = 

RHO( I 

P) 




I P 3 

= 

I PF ( 

IM— 1 , IT ) 




IP4 


IPF ( 

IM+l, IT) 




R ( 1 ) 

= 

RHO ( 

IP-1 ) 




IF ( 

IT 

.EQ. 

2) R ( 1 ) = 

RBV 

( IM 

, 1) 

R ( 2 ) 

= 

RHO ( 

IP-t-l) 




IF ( 

IT 

.EQ. 

I T M AX— 1 ) 

R( 2 

) = 

R BV ( I M , 2 ) 

R t 3 ) 

= 

RHO ( 

IP3) 




IF ( 

IM 

.EQ. 

2) R ( 3 ) = 

RBV 

(IT 

,3) 

R ( 4 ) 

= 

RHO ( 

IP4) 




IF ( 

IM 

.EQ. 

MMMM1 ) R 

( 4 ) = 

RBV ( IT, 4) 

BZ= 

BE 

( IM) 





B ( 3 ) 

= 

BE ( I 

M-I) 




B ( 4 ) 

= 

BE ( I 

M+l ) 





RETURN 

FNO 


SUBROUTINE AAK(IM,IP) 

AAK CALCULATES FINITE DIFFERENCE COEFFICIENTS, A, AND CONSTANT, K, 
AT A SINGLE MESH POINT 

COMMON / AUKRHO/ A { 2000 , 4 ) , U { 2000 ) , K ( 2000 ) , RHO ( 2 000 ) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM, MB I I M 1 , MB I IP 1 , MB00M1 , MBOOP 1 , MMMM1 , 
1HM 1,HM2,HM3,HT.DTLR, DMLR, PITCH. CP, EXPON, TV, W.CPTIP, TGROG, TBI.TBO, 
2LAMBDA,TWL,IT0R» IT MAX, NIP, IM S( 4 ) , B V ( 4 ) , MV ( 100 J , I V ( 10 1 ) , 

3UBVI 100,4) , RWB V ( 100, 4), I TV (100, 6), TV( 100 , 4 ) , DTUMV ( 100,4) , 

4BETAV( 100,4) ,MH( 100, 4 ) , DTDMH ( 1 00 , 4 ) , BET AH( 100,4) , RMH (100, 4) , 
5BEH( 100,4) ,RM{ 100) , BE (100) , D8DM ( 100 ) , SAL ( 100 ) , AAA ( 100) 

COMMON /HRBAAK/ H ( 4) , R ( 4 ) , B ( 4 ) , KAK ( 4 » , KA ( 4 ) , IH ( 4 ) , RZ , BZ 
INTEGER BLOAT, AANOK, ERSOR, STRFN, SURVL, AATEMP, SURF, 

IF IRST, UPPER, SI, ST, SRW 

REAL K, KAK, LAMBDA, LMAX.MH, MLE, MR , MSP , MV, MV I M 1 
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A12= 2./H( 1)/H(2) 

A34= 2./H{ 3)/H(4) 

AZ= A12+A34 

812= ( R (2 ) — R 1 1 J ) /RZZ t H( 1 )+H( 2) ) 

834= !B(4)*R(4)-B(3)*R(3) ) /BZ/RZ/! H(3)+H(4) ) -S AL ( I M ) /RM ( I M ) 
AlIP.l) = ( 2 • / H ( 1)+B12)/AZ/(H(1)+H(2)) 

A ( I P , 2 ) = A12/AZ-A!IP,1) 

A ( I P»3 ) = !2./H!3)+B34)/AZ/(H(3)+H(4) ) 

A ( I P » 4 ) = A34Z AZ-A ( I P » 3 ) 

K t IP) = -TWW*BZ*RZ*SAL( IM) /AZ 

RETURN 

END 


SUBROUTINE BDR Y 12 ( I , I M , I T , SURF ) 

BDRY12 CORRECTS VALUES COMPUTED BY HRB WHEN A VERTICAL MESH LINE 
INTERSECTS A BLADE 

COMMON /CALCON/ MB 1 1 , MBOO, MMM, MB I I Ml ,MB I I P 1 , MbOOMl , MROOP l , MMMM1 , 
1HM1 ,HM2. HM3.HT f DTLR, DMLR ,P I TCH , CP , EXPON , TW W, CP T I P , TGROG, TBI »TBO, 
2LAM3DA, TWL, I TOR , I TMAX , N I P. I MS ( 4 1 , B V ( 4 ) , MV ( 100 ) , I V ( 1 0 1 ) , 

3UBV (100,4) ,RWBV( 100,4 ) , I TV ! 1 00 , 6 ) , TV ( 100 , 4 ) , DTDMV! 100, 4) , 

4BETAV! 100 ,4) ,MH( 100, 4 ) ,D TDMH! 100, 4) , BET AH ( 100, 4 I ,RMH( 100,4) , 
5BEH( 100,4) ,RMC 100) ,BE ( 100) , DBDMI 100) ,SAL( 100) , AAA ( LOO) 

COMMON /RHOS/ RHOHBI 1 00, 4) , RHOVB ( 100 , 4 ) , RBV ( 100 , 4) 

COMMON / HR BAAK / HI 4 ) , R ( 4 ) , B ( 4 ) , KAK ( 4 ) » KA ( 4 ) , IH ( 4 ) , RZ , BZ 
INTEGER BLOAT, AANDK, ERSOR, STRFN, SUR VL , A ATE MP , SURF , 

1FIRST, UPPER, SL, ST, SRW 

REAL K • KAK , LAMBDA, LMAX,MH, ML F« MR ,MSP,MV»MVIM1 
H ( I ) = ABS ( FLOAT! IT-I TOR )*HT-TVl IM, SURF ) )*RM( IM) 

R { I )= RHOVB! I M, SURF) 

KAK(I)=BV!SURF) 

K A ( I ) = 1 

RETURN 

END 


SUBROUTINE BDRY34 ( I , I M , SURF ) 

BDRY34 CORRECTS VALUES COMPUTED BY HRB WHEN A HORIZONTAL MESH LINE 
INTERSECTS A BLADE 

COMMON /CALCON/ MBII ,MBOO,MMM,Mt)IIMI,MBI IPl,MbOOMl,MBOOPl»MMMMl» 
1HM 1»HM2?HM3,HT »DTLR» DMLR , PITCH, CP, EXPON, TWW»CPTIP» TGROG, TBI ,TBO, 
2LAMBDA,TWL» IT UR, IT MAX, NIP, IM S I 4 ) , B V { 4 ) , M V! 100 ) , I V ( 101 ) , 

3UBV(100,4),RWBV!100,4),ITV!100,6),TV(100,4),DTDMV(100,4), 

4BET AVI 100,4) , MH ( 100,4 ) ,DTDMHl 100,4) ,BETAH( 100,4) , RMH (100,4) , 
5BEH! 100 ,4) ,RM( 100) ,BE ( 100) , DBDM ( 100 ) , SAL ( 100 ) , AAA ( 100 ) 

COMMON /RHOS/ RHOHBI 100, 4) , RHOVB ( 100,4) , RBV ! 100, 4) 

COMMON /HRBAAK/ H ( 4) , R ( 4 ), B l 4 ), KAK ( 4 ), K A ( 4 ), I H ( 4 ), R Z , BZ 
INTEGER BLOAT, AANDK, ERSOR, STRFN, SUR VL , A A TEMP ,S URF , 

1 FIRST, UPPER, SI, ST, SRW 

REAL K , KAK,LAMBDA,LMAX,MH, MLE, MR , MSP , MV , MV I M 1 


56 



noon noon 


IHISURF)=IHISURF)+1 
I HS= IH ( SURF ) 

H ( I ) = ABS ( MV ( IM)-MH< IHS.SURF) ) 
RII )=RHDHB( IHS.SURF) 

BCE ) =BEHt IHS.SURF) 

KAK I I ) = BV (SURF) 

KA ( I ) = 1 

RETURN 

END 


SUBROUTINE SOR 

SOR SOLVES THE SET OF SIMULTANEOUS EQUATIONS FOR THE STREAM FUNCTION 

USING THE METHOD OF SUCCESSIVE OVER-RELAXATION 

COMMON / AUKRHO/ A ( 2000 , 4 ) , U { 2000 ) , K ( 2000 ) , RHO ( 2000 ) 

COMMON /INP/ GAM, AR,T IP, RHO IP, WTFL.WTFL SP.OMEGA.ORF.BE TA I , BE TAO, 
1N0BL.MBI , MB0,MBI2,MB02,MM»NBRI ,NBL,NRSP, MBDYF ,MBDYL»ITF, ITL, 
2BLDAT, AANDK,ERSOR,STRFN, INTVL, SURVL, MAGFAC , 
3MR(50)»RMSP150)»BESP150) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM, MB 1 1 MI , MB I I P 1 , MB00M1 , MBOOP l , MMMM 1 , 
LHMI ,HM2 ,HM3,HT ,DTLR, D MLR, PIT CH ,CP,EXPON,TWW,CPTIP, TGROG, TQI ,TBO, 
2LAMBDA.TWL, I TOR, ITMAX.NI P, IM S I 4 ) , BV I 4 ) , MV I 100) , I V( 101) , 

3UBVI 100,4) .RWBVI 100, 4 ) , I TV ( 1 00, 6 ) , T V ( 100 , 4 ) , DTDMV { 1 00 , 4 ) , 

4BETAVI 100,4) ,MHI 100, 4 ) , D TDMH ( 100 , 4 ) , BET AH ( 100,4) ,RMH( 100,4) , 

5BEH( 100,4 ) ,RM( 100) ,BE I 100) , DBDM ( 100 ) , SAL ( 100 ) , AAAI 100) 

INTEGER BLOAT, AANDK, ERSOR, STRFN , SURVL , AA TEMP , SURF , 

IF I R ST, UPPER, SI, ST, SR W 

REAL K , KAK , LAMBDA , LMAX,MH»MLE»MR,MSP»MV,MVIMl 
AATEMP = AANDK 
IF I0RF.GE.2.) ORF=0. 

IFIORF.GT.l. ) GO TO 20 
ORF = 1. 

ORFOPT = 2. 

10 ORFTEM = ORFOPT 
LMAX = 0. 

20 IFI AATEMP. GT.O) WRITE ( 6, 1010 ) 

ERROR = 0. 

SOLVE MATRIX EQUATION BY SOR, OR CALCULATE OPTIMUM OVERRELAXATION 

FACTOR 

DO SO I M=2 » MMMM1 
IFIAATEMP.GT.O) WRITE ( 6, 1020 ) IM 
IPU = IVI IM) 

I PL= I V ( IM+ 1 ) — 1 
IT= l 

DO 50 I P= I PU , I PL 
IFI IPU.GT.IPL) GO TO 50 
I T= 1T+1 

IFI IT.GT. I TV l IM,4) .AND. IT.LE.ITVI IM,3) ) 

1MAX0 1 1 TV I IM,4) , 1 ) — 1 
IFI IT.GT. I TV I IM,2) .AND. IT.LE. I TV I IM, 1 ) ) 

1MAX0I ITVI IM,2) , 1 )-l 


IT=IT+ITVIIM»3)— 
IT=IT+ITVl I M, 1 ) - 
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IP1= IP-1 
I P2= IP+1 
I P3= I PF ( I M— 1 , IT ) 

I P4— I PF ( I M+ 1 , I T ) 

I F ( IM.EQ.2) I P 3 = 0 
I F { IM.EQ.MMMM1 ) IP4 = 0 
IFlORF.GT.l. ) GO TO 30 
C CALCULATE NEW tSTIMATE FOR LM AX 

UNEW = A( IP, 1)*U( IP1 )+A{ IP,2)*U( IP2 )+A ( I P, 3 > *U ( I P3 ) +A < IP, 4) *U( IP4) 
IF (UNEW.LT.i.E-25) U(IP) = 0. 

IF (UUP) . FQ.O. ) GO TO 40 
RATIO = UNFW/UUP) 

LMAX= AMAX1 (RATIO, LMAX) 

UUP) = UNFw 
GO TO 40 

C CALCULATE NEW ESTIMATE FOR STREAM FUNCTION BY SOR 

30 CHANGE = ORE* ( K (I P ) -U ( I P ) + A ( I P , L ) *U (I PI ) +A ( I P , 2 ) *U (I P2 ) + A (I P , 3 ) * 
IU( I P3) +A ( I P , 4 ) *U ( IP4) ) 

ERRQR= AMAX1 ( ERROR, ABS(CHANGE) ) 

UUP) = U ( I P ) +CHANGE 
40 I F ( AATEMP.LE.O) GO TO 50 

WRITE (6,1030) IT, IP, IP1 ,IP2, IP3, IP4, ( A( IP , I ) , 1=1,4) ,K( IP) 

50 CONTINUE 
AA TEMP = 0 

IF(ORF.GT. 1- ) GO TO 60 
ORFOPT = 2./(l.+SQRT( ASS(l.-LMAX) ) ) 

WRITE (6,1040) ORFOPT 

IF (ORFTEM-URFOPT.GT. . 00001. OR. URFOPT.GT. 1.999) GO TO 10 
WRIIE (6,1000) 

ORF = ORFOPT 
GO TO 20 

60 IF(FRSOR.GT.O) WR I TE ( 6 , 1 050 ) ERROR 
I F ( ERROR. GT. .000001) GO TO 20 
IF( STRFN.LE.O) RETURN 

PRINT STREAM FUNCTION VALUES FOR THIS ITERATION 

WRITE (6,1060) 

I PL = 0 

MBOT = M I NO ( Mb I I M 1 , MB CO ) 

IF (MB0T.LT.2) GO TO 80 
DO 70 I M=2 » MBOT 
CALL WRITU( IM,5»2, IPL ) 

70 CALL WRITUt IM, 1,6, IPL ) 

80 MBI T = MAXO ( 2 , MB 1 1 ) 

MBOT = MINO ( MMMM1 , MBOO ) 

IF (MBIT. GT. MBOT) GO TO 100 
DO 90 I M=MBI T , MBOT 
CALL WR I TU ( I M, 5 , 4 , IPL ) 

CALL WRITU( IM.3,2, IPL ) 

90 CALL WRITUUM, 1,6, IPL ) 

GO TO 120 

100 MBIT = MAXO ( 2 , MB00P1 ) 

MBOT = MI NO ( MB 1 1 Ml »M MM Ml ) 

IF (MBIT. GT. MBOT) GO TO 120 
DO 110 I M=MB IT, MBOT 
110 CALL WRITUUM, 5,6, IPL ) 

120 MBIT = MAXOtMBII .MBOOPl) 
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IF IMBIT.GT.MMMM1 ) RETURN 
00 130 I M=MB IT, MMMM1 
CALL WRITUUM, 5,4, IPL ) 

130 CALL WRITUl IM,3,6» IPL ) 

RETURN 

1000 FORMAT ( 1H1 ) 

1010 FORMAT ( 8 2 HI IT IP IP1 IP2 IP3 IP4 All) A(2) 

1 A ( 3 ) A ( 4 ) K) 

1020 FORMAT ( 5H IM =,I4) 

1030 FORMAT ( IX, 14,516, 5F10.5) 

1040 FORMAT l 24H ESTIMATED OPTIMUM ORF =,F9.6) 

1050 FORMAT ( 8H ERROR =,Fll.8J 

1060 FORMAT ( 1H 1 , 10X»22HSTREAM FUNCTION VALUES) 

FND 


SUBROUTINE WRITUl I M, U PPE R, LOWER , IPL ) 

COMMON / AUKRHO/ A ( 2000 , 4 ) f U ( 200U ) , K l 2000 ) , RHO ( 2000 ) 

COMMON /CALCON/ MBII , MBOO, MMM , MB I I M 1 , Mb I I P l , MBOOM 1 , MBOOP 1 , MMMM1 , 
IHM1,HM2,HM3,HT,DTLR, DMLR ,P I TCH , CP , E XPON , TWW , CP TI P , TGROG, TB I , TBO, 
2LAMBDA , TWL » l TOR, ITMAX,NIP,IMS(4),BV(4),MV{100),IV(101), 

3UBV ( 100,4) , RWBV l 100, 4 ) , I TV( 100,6) ,TVl 100,4) ,DTUMV( 100,4) , 

4BETAVI 100,4) ,MH( 100,4) , DTDMH ( 100 , 4 ) , BE T AH ( 1 00 , 4 ) , R MH ( 1 00 , 4 ) , 
5BEHI 100,4 ) ,RM( 100) ,BE ( 100) , DBDM ( 100 ) , S AL { 1 00 ) , AA A { 1 00 ) 

INTEGER BLOAT, AANDK, ERSOR, STRFN , SUR VL , A ATE MP , S URF , 

1FIRST, UPPER, S1»ST»SRW 
I TVU = MAXOl ITVl IM, UPPER ), 2) 

I T VL = MINOlITVl IM, LOWER), ITMAX-i) 

I F ( ITVU.GT. ITVL ) RETURN 

IPU = IPL+1 

IPL = IPU+ITVL-ITVU 

WRITE16, 1000) I M , I TVU 

WR I TEC 6 , 10 10 ) (U( IP) , IP= IPU, IPL) 

RETURN 

1000 FORMAT ( 5H IM = , I 3 » 10X , 5H IT 1 =,I3) 

1010 FORMAT 12X,10F13.8) 

END 


SUBROUTINE SLAX 

SLAX CALLS SUBROUTINES TO CALCULATE RHO*W-SUB-M THROUGHOUT THE REGION 
AND ON THE BLADE SURFACES, AND TO CALCULATE AND PLOT THE 
STREAMLINE LOCATIONS 

COMMON /AUKRHO/ A l 20 00 , 4 ) , U ( 2000 ) , K ( 2000 ) , RHO l 2000 ) 

COMMON /INP/ GAM, AR, T I P , RHO I P , WTFL , WTFLSP , OMEGA , ORF , BETA I , BET AO, 

1 NOBL , MB I , M80,MBI2»MB02*MM, NB BI , NBL ,NRSP , MBDYF , MBDYL , I TF , I TL , 

2 BLOAT, AANDK, ERSOR, STRFN, INTVL , SURVL , M AGF AC , 

3 MR (50) ,RMSP(50) ,BESP( 50) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM, MB I I M l , MB I I P 1 , MBU0M1 , MB00P1 , MMMM1 , 
1HM 1, HM2, HM3, HT, DTLR, DMLR, PI TCH, CP, EXPON, TWW, CPTI P, TGROG, T8I, TBO, 

2 LAMBDA , TWL , I TOR , IT MAX, NIP, IMS{4),BV14),MV(100),IV1101), 

3UBVI 100,4 ) , RWB V ( 100, 4 ) , I TV l 1 00 , 6 ) , TV l 100 , 4) ,DTDMV ( 100,4) , 
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4BETAVC 100 ,4) ,MH( 100, 4 ) , DTDMH ( 100 , 4 ) , 8ETAH { 100,4) ,RMH{100,4) , 
5BEH( 100,4 ) ,KM( 100) .BE ( 100) ,D8DM( 100 ) , S AL ( 1 00 ) , AAA { 100) 
DIMENSION W (2000 ) ,RWM( 2000 ) , BETA ( 2000 ) , WMB (100,4) , WTB( 100,4) , 
lXDOWNt 8 00 ) , YACROSt 800 ) , T SL (400 ) , TSP( 100) ,USP ( 100 ) , DUDT (100) 
EQUIVALENCE ( A { 1 , 1 ) , W ( 1 ) ) , ( A ( 1 , 2 ) , RWM ( 1 ) ) , ( A < 1 , 3 ) , BET A ( 1 ) ) , 

1 (A( 1,4) , WMR ( 1 ) ) , (A(401,4),WTB{ 1) ) , ( A ( 80 1 , 4 ) , XOGWN ( 1 ) ) , 

2(K(1> , YACROS (1) ) , (K( 801) ,TSP( l) ) , 

3 ( K ( 901 ) ,USP{1) ), (KC1001) ,DUDT{ 1) ) 

INTEGER BLOAT, AANDK, ERSOR, STRFN, SURVL ,AA TEMP, SURF, 
IFIRST,UPPER,S1,ST,SRW 

REAL K, KAK, LAMBDA, LMAX.MH.MLE, MR, MSP, MV, MVIM1 

CALL SLAVBB THROUGHOUT REGION 

DO 10 1=1,4 
DO 10 I M= 1,1 00 
WMB ( I M , I ) = 0. 

10 WTB( IM, I ) = 0. 

MBOT = MI NO ( MB I I Ml , MBOO ) 

IF (MB0T.LT.2) GO TO 30 
DO 20 I M= 2 » MBOT 
CALL SLAVBB( IM.5,2) 

20 CALL SLAVBB( [M, 1,6) 

30 MBIT = MAXO ( 2 , MB I I ) 

MBOT = M I NO ( MMMM 1 , MBOO ) 

IF (MBIT. GT . MBOT ) GO TO 50 
DO 40 I M=M8 IT, MBOT 
CALL SLAVBBt IM.5,4) 

CALL SLAVBB ( IM, 3,2) 

40 CALL SLAVBBt IM, 1,6) 

GO TO 70 

50 MBIT = MAX0(2,MB00PL ) 

MBOT = MIN0(MBIIM1,MMMM1 ) 

IF (MBIT.GT.MBOT) GO TO 70 
DO 60 I M= MB IT, MBOT 
60 CALL SLAVBBt IM, 5, 6) 

70 MBIT = MAX01MBI I.MB00P1) 

IF (MB I T. GT. MMMM 1 ) RETURN 
DO 80 I M= MB I T , MMMM 1 
CALL SLAVBBt IM, 5,4) 

80 CALL SLAVBBt IM,3,6) 

RETURN 

END 
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SUBROUTINE SLAVBB ( IM , UPP ER , LOWER ) 

SLAVBB CALCULATES RHO*W-SUB-M ALONG VERTICAL MESH LINES 
COMMON SRW » ITER, IEND, L ER (2 ) , NER ID 

COMMON / AUKRHO/ A ( 2000 , 4 ) , U ( 2000 ) , K ( 2000 ) , RHO ( 2000 ) 

COMMON /INP/ GAM, AR, TIP, RHOI P,WTFL,WTFLSP, OMEGA, ORF, BE TAI .BETAO, 
1N0BL,MBI ,MB0,MBI2,MB02,MM,NBBI ,NBL , NRSP, MBOYF , MBOYL , I TF , I TL , 

2BLD AT, AANDK ,ERSOR,STRFN, INTVL , SURVL , MAGFAC , 

3MR ( SO) ,RMSP(50) ,BESP( 50) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM , MB 1 1 Ml , MB 1 1 P l , MB00M1 , MBOOP 1 , MMMMI , 
1HM1 , HM2 »HM3»HT ,DTLR, DMLR, PITCH, CP, EXPON,TWW,CPTIP, TGROG, TBI , TBO, 
2 LAMBDA, TWL, ITOR, IT MAX, NIP, IMS ( 4 1 , BV ( 4 ) , MV ( 100 ) , IV! 101) , 

3UBVI 100,4) ,RWBV< 100, 4 ) , I TV ( 1 00, 6 ) , TV ( 100 , 4 ) , DTDMV ( 100,4) , 

4BETAVI 100,4) ,MH( 100, 4 ) , DTDMH ( 100, 4 ) , BETAH! 100,4) ,RMH( 100,4) , 

5BEH { 100,4 ) ,RM{ 100) , BE ( 100) ,OBDM( 100 ) , SAL ( 100 ) , AAA ( 100) 

COMMON /RHOS/ RHOHBI 1 00 , 4) , RHOVB ( 100 ,4 ) , RB V ( 100, 4 ) 

DIMENSION W( 2000) ,RWM (2000) , BETA (2 000) , WMBI 100,4) ,WTB( 100,4) , 
1XD0WN! 800 ) . YACROSI800 ) , T SL ( 400 ) , TSP ( 1 00 ) , USP ( 100 ) , DUDT ( 100 ) 
EQUIVALENCE !A(1,1),W(1) ) , ( A ( 1 , 2 ) , RWM ( 1 ) ) , ( A < 1 ,3 ) , BETA ( 1 ) ) , 

1( A ( 1,4) ,WMR( 1) ), (A(401*4)*WTB( 1) ), { A ( 801 , 4 ) , XDQWN ( 1 ) ) , 

2 ( K ( 1 ) , YACROSI 1) ) , (KI801) , TSP ( 1 ) ), 

3 ( K ( 901 ) ,USP(1) ), (K(lOOl) ,DUDT( 1) ) 

INTEGER BLDAT, AANDK, ERSOR, STRFN , SURVL , AATEMP , SURF , 

IF IRST, UPPER, SI, ST, SRW 

REAL K, KAK, LAMBDA, LM AX, MH, ML E , MR , MSP , MV , MV I M 1 
I T VU= MAXOt 1TVI IM, UPPER) ,2) 

I TVL= MIN0( ITV( I M, LOWER ) ,ITM AX-1 ) 

NSP= ITVL-ITVU+3 
IFINSP.LT. 3) RETURN 
TSP ( 1 ) = FLOAT! 1-IT0R)*HT 

IF! ITV! IM, UPPERJ.LT. 2. OR. UPPER. EQ. 5) GO TO 10 
IF !TV!IM, UPPER). LT. TSP(l) ) GO TO 10 
TSP l 1) = TV!IM, UPPER) 

USP 1 1 ) = BV (UPPER ) 

GO TO 20 

10 USP ! 1 ) = UBV ( I M , 1 ) 

20 TSP(NSP)= FLOAT! ITMAX-ITOR )*HT 

IF! ITV! IM, LOWER) .GE. I TMAX. OR . LOWER . EQ .6 ) GO TO 30 
IF (TV! IM, LOWER). GE.TSP(NSP) ) GO TO 30 
TSP ( NSP ) = TV(IM, LOWER) 

USP ( NSP ) = BV ( LOWER ) 

GO TO 40 

30 USP ( NS P ) = UBV ( I M , 2 ) 

40 NSPM1= NSP— 1 
I T= 2 

IP= IPF ( I M , I TVU ) 

IPU= IP 

50 IF! IT.GT.NSPM1) GO TO 60 

TSP ! IT ) = FLOAT! IT-2+1 TVU-I TOR) *HT 
USP I IT ) = U( IP) 

I T= IT+ 1 
IP= IP+1 
GO TO 50 

CALCULATE RHO*W-SUB-M IN THE REGION, ANO RHO*W AT VERTICAL 
MESH LINE INTERSECTIONS ON THE BLADE SURFACES, OR RHO 
ON THE HORIZONTAL BOUNDARIES 

60 CALL SPLINE(TSP,USP, NSP, DUDT, AAA) 
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[ PL= IP-1 
I T= 2 
IP= IPU 

70 I F ( IP.GT. IPL) GO TO 80 

RWM { I P ) = DUDTl I T ) *WT FL/B E ( IMl/RMl I M ) 

IP= IP+1 
I T= IT+l 
GO TO 70 

80 IF( ITV( IM, UPPERI.LT. 2. OR. UPPER. £Q. 5) GO TO 90 
IF ( TV ( IM»UPPER) .LT. FLOAT! 1- I TOR ) *HT ) GO TO 90 
; UPPER BLADE SURFACE 

WMB ( IM,UPPER)= DUDTl 1 ) *W TFL/ BE ( I M ) /RM ( I M ) 

RMDIU2 = (RM(IM)*DTDMV( IM, UPPER) )**2 
I F ( RMDTU2.GT. 10000.) WMB ( I M, UPPER ) =0. 

WMB ( IM, UPPER) = A8S( WMB( IM, UPPER) )*SQRT { 1. + RMDTU2 ) 

GO TO 100 
C LOWER BOUNDARY 

90 RWMbV = DUDr (1)*WTFL/BE( IM)/RM( IM) 

RW= SORT! RWBV ( IM, 1 ) * * 2+R WMBV ** 2 ) 

TWLMR= 2. *OMEGA*LAMBDA-( 0MEG4*RM( IM) )**2 
LERI 1 ) = l 

C OENSTY CALL NO. 1 

CALL DENSTY l RW , RBV ( I M, 1) , ANS , TWLMR , CPT I P , E XPON , RHO I P , GAM , AR , TI P ) 
100 I F ( I TV ( IM, LOWER). GE. ITMAX. OR. LOWER. EQ. 6) GO TO 110 
IF ( TV I IM, LOWER) .GE. FLOAT! I TMAX- I TOR ) *HT ) GO TO 110 
C LOWER BLADE SURFACE 

WMB I IM, LOWER ) = UUDTI NSP ) *WTFL/BE I IM) /RMI IM) 

RMDTL2 = (RMI IM)*DTDMV( IM, LOWER) )**2 
IF (RMDTL2.GT. 10000.) WMB 1 1 M, LOWER ) =0. 

WMB I IM, LOWER ) = ABSIWMRI IM, LOWER) ) *SQRT ( 1.+RMDTL2) 

RETURN 

C UPPER BOUNDARY 

110 RWMI.W = DUDT(NSP)*WTFL/BE( IMl/RMI IM) 

RW= SORT (RWBV l I M , 2 )* *2+R WMBV ** 2 ) 

TWLMR= 2. * OMEGA* L AMB DA- ( OMEGA* RM( IM) )**2 
LER(l) = 2 

C DENSTY CALL NO. 2 

CALL DENSTY(RW,RBV(IM,2) , ANS , TWLMR , C PT I P , E XPON , RHO I P , GAM , AR , T I P ) 

RETURN 

END 


SUBROUTINE TAnG 
C 

C TANG CALCULATES RHO*W-SUH— TEET A AND THEN RHO*W THROUGHOUT THE REGION 
C AND UN THF BLADE SURFACES, AND CALCULATES THE VELUCITY ANGLE, BETA, 

C THROUGHOUT THE RtGION 
C 

COMMON SRW, ITER, I LND, LER (2 ) , NER( 1 ) 

COMMON /AUKRHU/ A( 2000,4 ),U( 2000) ,K( 2000) ,RHO( 2000) 

COMMON /INP/ GAM, AR,T IP, RHO IP, W TFL, W TFL SP, OMEGA, OK F, BETA I , BE TAD, 

1 NOBL ,MB I , M p 0 , M b I 2 , MB 0 2 , M M , N B B I , NRL , NRSP , MBDYF , MBDYL , ITF, ITL, 
2BLDAT , AANDK , EKSOR, ST REN, INTVL, SURVL , MAGE AC , 

3MR ( >0) , R M S P ( SO) ,BESP ( 50) 
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COMMON /CALCON/ Mbll , MBOO, MMM , MB I I M 1 , MB I I P 1 , MbOOMl , MBOOP 1 , MMMM1 , 
1HMI ,HM2,HK3,HT,OTLR,OMLR,P ITCH ,CP , E XPON , TWW , CPT I P , TGROC , TBI ,TBO, 
2LAMBDA, TWL, l TOR, UMAX, NIP, IMS(4),BV(4),MV(100),IV(101), 

3UBV ( i00»4 ) , RWbV ( 100, 4) , I TV ( 1 00 , b ) , TV ( 1 00 , 4 ) , OTOM V ( 100,4) , 

43ET4VI 100,4) ,MH( 10 0,4) ,DTDMH( 100,4 ) , BET AH ( 1 00 , 4 ) , R MH ( 1 00 , 4 ) , 

5BEH ( 100,4) ,RM( 100) , BE ( 10 0) ,D3DM{ 100) ,SAL (100) , AAAI 100) 

COMMON /RHOS/ RH0H3 ( 1 00 , 4 ) , RHO Vt5 ( 100 , 4 ) , RB V ( 1 00 , 4 ) 

DIMENSION SPM{ 100) ,USP( 100) , DUQMC 100) 

DIMENSION W (2000) , RWM ( 2000 ) , RE TA ( 2000 ) , WMB ( 100 ,4 ) , WTB ( 100 ,4 ) , 
lXDO'wNM 800 ) , Y ACROS ( 00 0) 

EQUIVALENCE ( A ( 1 , 1 ) , W { 1 ) ) , ( A ( 1 , 2 ) , RWM ( 1 ) ) , ( A ( 1 , 3 ) , BE TA ( 1 ) ) , 

1 ( A( 1 ,4) , WMB ( 1 ) ) , ( A (40 1,4 ),WTB( 1 ) ) , ( A( 801,4) ,XDGWNl 1 ) ) , 

2 ( K ( 1 ) , YACROS ( 1 ) ) 

INTEGER BLOAT, AANDK, ERSOR, STRFN, SURVL , AATEMP , SURF, 
1FIRST,UPPER,S1,ST,SRW 

REAL K,KAK, LAMBDA, LMAX,MH,MLE, MR, MSP, MV, MVIM1 
EXTERNAL BL 1 , BL2 , BL3 , BL4 

PERFORM CALCULATIONS ALONG ONE HORIZONTAL LINE AT A TIME 
I T= 2 

10 IF ( IT.EQ. UMAX) RETURN 

CN GIVEN HORIZONTAL MESH LINE, FIND FIRST POINT IN THE REGION 

IF ( MB II. LE. 1. AND. (IT. LE. ITV (1,4). OR. (IT. GE.ITVtl, 3). AND. 

1IT.LE. ITV( 1,2) ) . OR . IT ,GE.ITV( 1, l) ) ) GO TO 50 
IF (MB I I ,GF. 2. AND. ( MBCO . LT . 0 .OR . I T . L E . IT V ( 1 , 2 ) -OR . 

UT.bE. ITV( 1,1) ) ) GO TO 50 
IM = 1 
20 IM= IM+1 

IF (IM.GE.MMM) GO TO 180 

00 30 SURF= 1,3,2 

IF ( IM.GT.MBCO. AND. SURF. EQ. 1 ) GO TO 30 
IF ( IM.LE.MbI I .AND.SURF.EQ.3) GO TO 30 

IF ( IT.GE. ITV ( IM.SURF ) .AND. IT.LT. ITV( IM— 1 , SURF ) ) GO TO 60 
30 CONTINUE 
SURF = 1 

IF ( IM. EQ. MBCoP 1. AND. IT . EQ . I TV IMBOO, 1 )- L.AMO. 

1 ITV (MBOO, 1 )-ITV ( MBOO, 2 ) . EQ.2 ) GO TO 60 
DO 40 SURF=2»4»2 

IF ( IM.GT . M BOO. AND. SURF. EQ.2 ) GU TO 40 
IF ( IM.LE.MBI I. AND. SURF. LQ. 4) GO TO 40 

IF ( IT. LE . IT V ( IM, SURF ) . AND. IT.GT. ITV ( I M- 1 , SURF ) ) GO TO 60 
40 CONTINUE 
GO TO 20 

FIRST POINT IS ON LEFT BOUNDARY 

50 SI = 0 
I M 1 = 1 

1 M= 2 

SPM ( 1 ) = M V ( 1 ) 

USP( I)= UBV ( IT, 3) 

GO TO 70 

FIRST POINT IS ON A BLADE SURFACE 


60 S 1= SURF 
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I M 1 = I M— 1 
IM2 = IM 

TH= FLOAT ( T T — I TOR ) *HT 
MVIMl = M V ( INI) 

IF (IM.EQ.MbIIPl.AND.(SURF.EQ.3.0R.SURF.EQ.4)) MVIM1= 

1 MV I M 1+ ( MV ( I M2 ) —MV I M 1 )/IOOO. 

LER ( 2 ) = 9 

C BLCD (VIA RGOT ) CALL MO. 9 

IF (Sl.EQ. I. AND. IMl.NE.MbOQ) CALL ROOT ( MVI Ml , MV { I M2 ) , TH , BL 1 , 
1DTLR » ANS t A Aa ) 

LER ( 2) = 10 

C 8LCD (VIA ROOT) CALL NO. 10 

IF (Sl.EQ. 2) CALL ROO T ( M VI M 1 , MV ( I M2 ) , TH , BL2 , D TLR , ANS , AAA ) 

LER ( 2 ) = 11 

C BLCD (VIA ROOT) CALL NO. 11 

IF (Sl.EQ. 3) CALL ROOT ( MVI M 1 ,MV l IM2 ) , TH , BL 3, DTLR , ANS , AAA ) 

LER ( 2) = 12 

C BLCD (VIA ROOT) CALL NO. 12 

IF (Sl.EQ. A) CALL ROOT (M VIM 1 ,MV ( I M2 ) » TH , BL4, DTLR , ANS, AAA ) 

IF (Sl.EQ. 1. AND. IM1. EG. MHOO) ANS=MV(MBOO) 

SPM ( I Ml ) = ANS 
US P ( I M 1 ) = BV(Sl) 

C 

C MOVE ALONG HORIZONTAL MESH LINE UNTIL END OF REGION IS REACHED 

C 

70 DO dO SURF= 1.3.2 

IF ( IM.GT.NBOU. AND. SURF. EQ.l ) Gu TO 80 
IF (IM.LE.Miil I. AND. SURF. EQ. 3) GU TO 80 
IF (ITVIIM-1, SURF). EQ. -10000) GU TO 80 

IF ( IT.LT. I FV( IM»SURF ) .AND. IT.GE. ITV( IM-1 , SURF ) ) GO TO 110 
80 CONTINUE 
SURF = 3 

IF ( IM.EQ.MBl I .AND. IT. EQ.l TV (MB I I, 3)- 1. AND. 

1ITVIMBI I ,3)-ITV(MBII , 4 ) . EQ. 2 ) GO TO 110 
DO 90 SURF=2»4»2 

IF { IM.GT.Mu.OO.AND.SURF.FQ.2) GO TO 90 
IF ( IM.LE.MBl I .AND. SURF. EQ.4) GO TO 90 

IF ( IT. GT. ITV( IM, SURF ) .AND. IT. Lt. ITV( IM-1, SURF ) ) GO TO 110 
90 CONTINUE 

SPM ( IM ) = MV( IM) 

IP= IPF(IM.IT) 

USP( I M ) = U( rp) 

IF (IM.EQ.MMM) GO TO 100 
I M= IM+1 
GO TO 70 
C 

C FINAL POINT IS ON RIGHT LOUNDARY 

C 

100 ST = 0 

IMT = MMM 

USP(IMT) = UBV I IT, 4) 

GO TO 120 
C 

C FINAL POINT IS ON A BLADE SURFACE 

C 

110 ST= SURF 
I MT= IM 
I MTM1= IMT— 1 
TH= FLOAT ( IT-ITOR) *HT 
MV I Ml = MVI1MTM1) 
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IF ( (IMTM1.EQ.MBII ). AND. ( S T . EQ . 3 . OR. ST. EQ. 4 ) . AND . ( I TV (MB I I ,3 >- 
1 I TV ( MB I I i 4 ) . EQ.-2 ) ) MVIM1 = M V I M 1+ ( M V { I M T ) -M V I M 1 ) / 1 000 . 

L ER ( 2 ) = 13 

C BLC 0 (VIA ROOT) CALL NC. 13 

IF (ST.EQ.l) CALL ROO T { MVI Ml f M V ( I MT ) , TH , BL 1 , OT LR , ANS , A AA ) 

LER { 2 ) = 14 

C BLCC (VIA ROOT) CALL NO. 14 

IF (ST.FQ.2) CALL ROOT ( M VI Ml , MV ( IMT ) , TH , BL 2 , DTLR , ANS , A AA ) 

LER ( 2 ) = 14 

C BLCO (VIA ROOT) CALL NO. 15 

IF (ST.EQ.3.AND.IMT.NE.M0II) C A lC'-'R-O 0 T ( MVI M 1 , M V l 1 M T ) , TH. BL3» 
10TLR.ANS, AAA) \ 

LER ( 2 ) = 14 

C BLCD (VIA ROGT) CALL NO. 16 

IF (ST.E0.4) CALL ROOT ( MVI Ml ,MV ( IMT ) . TH, BL4, OTLR , ANS , AAA ) 

IF (ST. EO. 3. AND. IMT.EC.MI3I I ) ANS=MV(MBII) 

SPM ( IMT ) = ANS 
US P ( I M T ) = 8 V ( S T ) 

CALCULATE RHL’*W — SUB- THETA AND THEN RHC*W AND BETA IN THE REGION 
120 NS P= IMT-IMi+l 

CALL SPLI NE ( SPM( I Ml) ,USP( I Ml ) , NSP, DUDM( I M 1 ) , AA A ( I M 1 ) ) 

F I R ST= 2 

IF (IM1.NE.1) FIRST=IM2 
LAST = MMMM 1 

IF (IMT.NE.MMM) LAST= IMTM1 
IF (FIRST. Gl .LAST) GO TO 140 
DO 130 I=FIRST , LAST 
RWT = -OUDM( I )*WTFL/BE( I ) 

IP = I PF ( I , I T ) 

W ( IP)- SORT (RWT»*2+RWM( I P) **2 ) 

130 (3 E T M I P ) = ATAN ( RWT/RWM ( IP) )*57. 295779 

CALCULATE RHO*W UN THE BLADE SURFACES, OR RHO ON VERTICAL BOUNDARIES 

140 IF (Sl.EQ.O) GO TO 150 

CALL SEARCH ( S P M ( I M 1 ) , S 1 , I H S ) 

ANS- -DUDM ( IM1)*WTFL/BEH(IHS,S1) 

WTB(IHS,S1)= A3SIANS) *SQRT( l.+l./(RMH( I H S , S 1 ) * DTDMH ( IHS,S1>)**2) 

GO TO 160 

150 RWT= -DOOM ( 1 ) *WTFL/8E ( 1) 

RW= SQRT ( RWT « * 2 + R W S V ( IT,3)**2) 

TWLMR= 2. *OMEGA*LAMGDA-( CMEG A*RM ( 1 ) ) **2 
LER(l) = 3 
C DENSTY CALL NO. 3 

CALL DLNSrY(RW,RBV(IT ,3) , AN S , THL MR , C PT I P , E XPON , RHO I P , GAM , AR , T I P ) 
160 IF (ST.EQ.O) GO TO 170 

CALL SEARCH(SPM( IMT) , ST, IHS) 

ANS= -DUDM ( IMT )*WTFL/BEH( IHS, ST) 

WTB( IHS , S T ) = ABS(ANS) *SQRT( l. + l./(RMH( I HS , ST ) * DTDMH ( IHS,ST))**2) 
GO TO 20 

170 RWT= -DUDM ( MMM ) *W TFL / B E ( MMM ) 

RW= SQRT(RtoT**2+RWBV( IT,4)**2) 

TWLMR= 2 . *OMEGA*L AMBDA- ( OMEGA* RM ( MMM ) )**2 
LER ( 1 ) = 4 

C DENSTY CALL NO. 4 

CALL DENSTY (RW.RBVl I T , 4 ) , ANS , T WLMR , C P T I P , E XPON , RHO I P , GAM , AR , T I P ) 
180 I T= IT+1 
GO TO 10 
END 
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SUBROUTINE SEARCH ( D I ST, SURF , I S ) 

C 

C SEARCH LOCATES THE POSITION OF A GIVEN VALUE OF M IN THE MH ARRAY 
C 

COMMON /CALCQN/ MBII , MBOO, MMM, MB 1 1 MI , MB I IP 1 , MBOOM1 , MBOOP1 , MMMM1 , 
IHMI , HM2, HM3, HT , DTLR, OMLR , PITCH, CP, EXPON , TWW, CPTI P , TGROG , TBI ,TBO, 
2LAMBDA,TWL»IT0R, IT MAX, NIP, IM S ( 4 1 , B V ( 4 ) , MV ( 100) , I V( 1011 , 

3UBV ( 100,4) .RWBVI 100, 4 ) , I TV ( 100, 6 ) , T V ( 100 ,4 ) ,DTDMV ( 100,4) , 
4BETAV(100,4) ,MH( 100, 4) ,DTDMH{ 100,4) ,BETAH( 100,4) ,RMH( 100,4) , 

5REH( 100,4) ,RM{ 100) ,BE ( 100) , DBOM { 100 ) , SAL { 1 00 ) , AAA ( 100 ) 

INTEGER BLOAT, AANDK, ERSOR, STRFN, SURVL , AATEMP , SURF, 

1FIRST, UPPER, Sl.ST.SRW 

REAL K,KAK,LAMB0A,LMAX,MH»MLE»MR,MSP,MV»MVIM1 
DO 10 1=1,100 

IF ( ABS I MH { I , SURF )— D 1ST) .GT . OMLR ) GO TO 10 
IS = I 
RETURN 
10 CONTINUE 

WRITE (6, 1000) 0 I ST, SURF 
STOP 

1000 FORMAT (38HL SEARCH CANNOT FIND M IN THE MH ARRAY/7H DIST =,G14.6, 
1 10X , 6HSURF =,G14.6) 

END 


SUBROUTINE VELOCY 
C 

C VELOCY CALLS SUBROUTINES TO CALCULATE DENSITIES AND VELOCITIES 
C THROUGHOUT THE REGION AND ON THE BLADE SURFACES, AND IT PLOTS 
C THE SURFACE VELOCITIES 
C 

COMMON / AUKRHQ/ A ( 2000 ,4 ) , U ( 2000 ) , K ( 2000 ) , RHO ( 2000 ) 

COMMON /INP/ GAM, AR, TIP, RHOI P » W TFL » WTFL SP, OMEGA, OR F ,BETAI , BET AO, 
1N0BL.MBI ,MBG,MBI2,MB02,MM,NBBI , NBL , NRSP , MBDYF , MBDYL , ITF, ITL, 
2BLDAT, AANDK, ERSOR, STRFN, INTVL , SURVL , MAGFAC , 

3MR(50) , RMSPC50) ,BESP(50) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM , MB 1 1 Ml , MB I I PI , MB00M1 , MB OOP 1 ,MMMM1 , 
1HM1,HM2,HM3,HT»DTLR» OMLR ,P I TCH , CP , EXPON , TWW, CP T I P , TGROG, TB I , TBO, 
2LAMBDA , TWL , I TOR , I TMAX , N I P, IMS(4),BV(4),MV(100),IV1101>, 

3UBV l 100,4) ,RWBV( 100,4) , ITV( 100, 6), TV ( 100,4) ,DTDMV( 100,4) , 

4BETAV( 100,4) ,MH( 100,4 ) , 0 TDMH ( 100, 4 ) , BET AH( 100, 4 ) , RMH ( 100, 4 ) , 
5BEH( 100,4) ,RM( 100) ,BE ( 100) , DBDM( 1 00 ) , SAL ( 1 00 ) , AAA ( 100) 

DIMENSION KKK ( 18 ) 

DIME NS ION W ( 2000 ) , RWM { 2000 ), BETA ( 2000 ), WMB ( 1 00 ,4 ), WTB ( 100,4 ) , 

1 XDOWNl 800 ) , Y ACROS ( 800 ) 

EQUIVALENCE (A(l,l),W(l)), ( A ( l , 2 ) , RWM ( 1 ) ) , ( A ( 1 ,3 ) , BETA ( 1 ) ) , 

1 (A ( 1,4) , WMB ( 1 ) ) , (A(401,4),WTB( 1) ) , ( A ( 80 1 ,4 ) , XDUWN ( 1 ) ) , 

2 I K ( 1 ) , YACROS ( 1 ) ) 

INTEGER BLOAT, AANDK, ERSOR, STRFN, SURVL , AATEMP , SURF , 

1FIRST, UPPER, S1,ST,SRW 

REAL K,KAK,LAMBOA,LMAX,MH, ML E, MR, MSP , MV , MV I M 1 

DATA KKK( 4)/ 1H*/ ,KKK( 6) / 1H0/ , KKK { 8 ) / 1H=/ , KKK ( 10)/1H{/, 

1KKK ( 12) /1H+/.KKKI 14) / 1HX/,KKK{ 16 ) / 1H$/ , KKK ( 1 8 ) / 1H ) / 

C 

C CALL VELBB AND VELSUR THROUGHOUT THE REGION 

C 
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MBOT = MI NO ( MB 1 I Ml » MBOO ) 

IF (MB0T.LT.2) GO TO 20 
00 10 IM=2,MBUT 
CALL VELBBI IM, 5,2) 

10 CALL VELBBI IM, 1,6) 

20 MBIT = MAXO ( 2 , MBI I ) 

MBOT = MIN0(MMMM1,MB00) 

IF (MBIT. GT. MBOT) GO TO 40 
00 30 i M=MB I T, MBOT 
CALL VELBBI IM, 5, 4) 

CALL VELBBI IM, 3,2) 

30 CALL VELBBI IM, 1,6) 

GO TO 60 

40 MBIT = MAXO ( 2 , MB00P1 ) 

MBOT = MINOIMBI IM1,MMMM1 ) 

IF (MBIT. GT. MBOT > GO TO 60 
DO 30 I M=MB IT, MBOT 
50 CALL VELBBI IM, 5,6) 

60 MBIT = MAXOIM0I I .MB00P1 ) 

IF (MBIT. GT . MMMM1 ) GO TO 80 
DO 70 I M=MB I T , MMMM1 
CALL VELBBI IM,5,4) 

70 CALL VELBBI IM,3,6) 

80 CALL VELSUR 

PREPARE INPUT ARRAYS FOR PLOT OF VELOCITIES 
NP2= 0 

C SURFACES 1 TO 4 - TANGENTIAL COMPONENTS 
DO 110 SURF = 1,4 
NP1= NP2 
IMSS= I MS I SURF ) 

IF (IMSS.LT.l) GO TO 100 
DO 90 I HS=1 , IMSS 

IFIWTBI IHS, SURF). EQ.O. ) GO TO 90 

IF (ABSIDTDMHI IHS , SURF ) *RMH( IHS, SURF) ).LT.. 57735) GO TO 90 
NP 1 = NPl+l 

YACROS I NP I ) = WTBI IHS, SURF) 

XDQWNINPl ) = MHIIHS.SURF) 

90 CONTINUE 

100 KKK(2*SURF+l )= NP1-NP2 
110 NP2= NP 1 

C SURFACES l AND 2 - MERIDIONAL COMPONENTS 
DO 140 SURF= 1,2 
NP 1= NP2 

MBOT = MI NO ( MBOOM 1 »MMMM1 ) 

IF I 2.GT .MBOT ) GO TO 130 
DO 120 I M=2 , MBOT 

IFIWMBI IM,SURF) .EQ.O. ) GO TO 120 

IF (ABSIDTOMVI IM,SURF)*RM( IM) ).GT. 1.7321) GO TO 120 
NP1= NP 1+ 1 

YACROS I NP 1 ) = WMB I I M, SURF ) 

XDOWNI NPl ) = MVIIM) 

120 CONTINUE 

130 KKK I 2* SURF +9 ) = NPI-NP2 
140 NP2= NPl 

C SURFACES 3 AND 4 - MERIDIONAL COMPONENTS 
DO 170 SURF=3»4 
NP1= NP2 

MBIT = MAXO I MB I I P 1 , 2 ) 


IF(MBIT.GT.MMMMl) GO TO 160 

DO 150 IM=M6IT,MMMM1 

I F { WMB ( IM, SURF). EQ. 0.1 GO TO 150 

IF ( ABS(DTDMV( IM,SURF )*RM( IMH.GT. 1.7321)- GO TO 150 
NP 1= MP 1+1 

YACROS ( NP 1 ) = WMB(IM,SURF) 

XDOWN( NP1 ) = MV(IM) 

150 CONTINUE 

160 KKK ( 2* SURF + 9 ) = NP1-NP2 
170 NP2= NP 1 
C 

C PLOT VELOCITIES 

C 


1000 

1010 


KKK ( l ) = 1 
KKK 12 ) = 8 
P= 5. 

WR I TE ( 6 , 1000 ) 

CALL PLOTMY( XDOWN , YACROS , KKK , P } 
WR I TE ( 6 , 1 0 10 ) 

RETURN 

FORM AT (2HPT , 50X , 24HBL ADE SURFACE 
FORMAT (2HPL,37X»63H VELOCITY (W) 
1(M) DOWN THE PAGE / 2HPL / 


VELOCITIES) 
VS. MERIDIONAL 


STREAMLINE DISTANCE 


22HP L , 50 X , 

50H+ - 

BLADE 

SURFACE 

1, 

BASED 

ON 

32HPL, 50X, 

50H* - 

BLADE 

SURFACE 

1, 

BASED 

ON 

42HPL.50X, 

50HX - 

BLADE 

SURFACE 

2, 

BASED 

ON 

52HPL ,50X, 

50H0 - 

BLADE 

SURFACE 

2, 

BASED 

ON 

62HPL »50X , 

50HS - 

BLADE 

SURFACE 

3, 

BASED 

ON 

72HPL.50X, 

50H= - 

BLADE 

SURFACE 

3, 

BASED 

ON 

82HPL ,50X , 

50H ) - 

BLADE 

SURFACE 

4, 

BASED 

ON 

92HPL.50X, 

5 OH ( - 

BLADE 

SURFACE 

4, 

BASED 

ON 


MERIDIONAL 

TANGENTIAL 

MERIDIONAL 

TANGENTIAL 

MERIDIONAL 

TANGENTIAL 

MERIDIONAL 

TANGENTIAL 


COMPONENT/ 

COMPONENT/ 

COMPONENT/ 

COMPONENT/ 

COMPONENT/ 

COMPONENT/ 

COMPONENT/ 

COMPONENT) 


END 


SUBROUTINE V ELBB ( I M, U PPE R, LOWER ) 

C 

C VEL CALCULATES DENSITIES AND VELOCITIES FROM THE PRODUCT OF 
C DENSITY TIMES VELOCITY 
C 

COMMON SRW, ITER, I END, LER (2 ) , NERC 1) 

COMMON / AUKRHO/ A ( 2000 , 4 ) , U ( 2000 ) , K ( 2000 ) , RHO { 2000 ) 

COMMON /INP/ GAM, AR, T I P , RHO I P , WTFL , WTFLSP , OMEGA, ORF , BETA I , BETAO, 
1N0BL,MBI»MB0,M6I2»MB02,MM,NBBI, NBL »NRSP»MBDYF, MBDYL »ITF» ITL, 
2BLDAT, AANDK.ERSOR, STRFN, INTVL , SURVL , MAGF AC , 

3MR150) , RMSP (50) ,BESP( 50) 

COMMON /CALCON/ MBII , MBOO, MMM , MBI IM1 , MB I I P 1 , MB00M1 , MB00P1 ,MMMM1 , 
1HM1,HM2»HM3,HT ,DTLR, DMLR ,P I TCH, CP, EXPON , TWW, CP TI P, TGROG, TBI , TBO, 
2LAMBDA»TWL, IT OR, IT MAX ,NIP, IMS(4),BV(4),MV( 100) ,IV( 101) , 
3UBVKOO,4) ,RWBV( 100,4) , I TV ( 100,6) , TV( 100, 4 ) , DTDM V { 1 00, 4) , 

4BETAVI 100,4) ,MH( 100, 4 ) , DTDMH { 100, 4 ) , BETAH { 100, 4 ) , RMH { 100,4 ) , 
5BEH( 100,4) ,RM( 100) ,BE ( 100) ,DBDM( 100) , SAL { 100 ) , AAA ( 100 ) 

COMMON /RHOS/ RHOHBI 1 00 , 4) , RHOVB I 100, 4 ) , RB V { 100 , 4) 

DIMENSION WWCRM( 100,4) ,WWCRT t 100,4) 

DIMENSION W ( 2000 ) ,RWM( 2000) , BETA ( 2000 ) , WMB ( 1 00 , 4 ) , WTB ( 100 , 4 ) , 
1XD0WNI 800), YACROS (800) 

EQUIVALENCE (A(1,1),W(1)), ( A ( 1 , 2 ) , RWM ( 1 ) ), ( A ( 1,3) , BETA ( 1 ) ) , 

1 (A (1,4) ,WMB(1) ) , (A(401,4),WTB{ 1) ), ( A ( 801 , 4) , XDOWN ( 1 ) ) , 
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2 (K( l ) , YACROS (1)1 

INTEGER BLOAT, AANDK, ERSOR, STRFN»SURVL,AA TEMP, SURF, 

1FIRST, UPPER, Sl.ST, SR W 

REAL K.KAK, LAMBDA, LMAX,MH,MLE, MR, MSP, MV, MVIM1 

C VELBB CALCULATES ALONG VERTICAL MESH LINES FROM BLADE TO BLADE 

C 

IF ( IM. NE . 2. OR. UPPER . NE . 5) GO TO 10 
IF{ INTVL.GT.O) WR ITE l 6, 1 000 ) 

RELER = 0. 

10 ITVU =MAXO( ITV( IM, UPPER) ,2) 

ITVL = MINO( ITV( IM, LOWER ), ITMAX-1 ) 

IPUP1 = IPFIIM.ITVU) 

IPLM1 = IPFdM.ITVLJ 

TWLMR= 2 . *OMEGA*L AMBD A- ( OMEGA* RM ( IM) )**2 
WCR= SORT ( T GROG*T I P* ( l.-TWLMR/CPTIP ) ) 

IF ( ITVL. LT. ITVU) GO TO 30 
C ALONG THE LINE BETWEEN BLADES 

DO 20 IP= I PUP 1,1 PLM 1 
LER(l) = 5 

C DENSTY CALL NO. 5 

CALL DENSTY (W( IP) ,RHO{ IP ), ANS , TWLMR.CPT I P, EXPUN, RHOI P ,GAM, AR.TIP) 
20 W(IP)= ANS 

IF (INTVL.LE.O) GO TO 30 

WRITE (6, 10 10) IM, (Wt I P ) , BET A ( I P ) , I P= IPUP1 , IPLM 1 ) 

C ON THE UPPER SURFACE, IF IT IS A BLADE 

30 IF (UPPER. EQ.6) GO TO 40 

IF lITV(IM,UPPER).LT.2) WMB ( IM , UPPER ) = 0. 

RHOP.= RHOVB( IM, UPPER ) 

LER(l) = 6 

C DENSTY CALL NO. 6 

CALL DENS TY( WMB ( IM.UPPER ) , RHOVB( IM , UPPER ), ANS , TWLMR ,CPTI P , EXPON, 
1RH0IP, GAM, AR, TIP) 

WMB ( IM , UPPER ) = ANS 

WWCRM( I M , UPPER )= WMB ( IM, UPPER > /WCR 

RELFR= AMAX1 (RELER, ABSl( RHOB— RHOVB l IM.UPPER) ) / RHOVB ( I M , UPPER ) ) ) 
C ON THE LOWER SURFACE, IF IT IS A BLADE 

40 IF (LOWER. EQ. 6) RETURN 

IF (ITVIIM, LOWER). GT. ITMAX-1) WMB ( IM, LOWER ) = 0. 

RH08= RHOVtH IM, LOWER) 

LER(l) = 7 

C DENSTY CALL NO. 7 

CALL DENSTY ( WMB ( IM, LOWER), RHOVB! I M, LOWER ) , ANS, TWLMR.CPT IP, EXPON, 
1RH0IP, GAM, AR.TIP) 

WMB( IM,LOWER)= ANS 

WWCRM! IM, LOWER ) = WMB ( I M , LOWER ) /WCR 

RELER= AMAX1(RELER»ABS( ( RHOB-RHOVB ( IM, LOWER) ) / RHOVB ( I M, LOWER )) ) 
RETURN 

VELSUR CALCULATES ALONG A BLADE SURFACE 

ENTRY VELSUR 
DO 60 SURF=1 ,4 
I MSS = IMS ( SURF ) 

IF (IMSS.EO.O) GO TO 60 
DO 50 IHS=1, IMSS 

TWLMR= 2. *OMEG A*L AMBDA- { QMEGA*RMH( IHS.SURF ) )**2 
WCR= SORT ( TGROG*T I P* (l.-TWLMR/CPTIP) ) 

RHOB= RHOHB( IHS.SURF ) 

LER(l) = 8 
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C OENSTY CALL NO. 8 

CALL DENSTY! WTB! IHS, SURF ) , RHOHB ( I H S, SURF ) , ANS, TWLMR.CPTIP.EXPON, 
1RHOIP.GAM, AR.TIP) 

WT8 ( IHS , SURF ) = ANS 

kWCRT ( IHS* SURF ) = WTB( I HS ,SUR F ) / WCR 
50 RELER= AMAX1(R£LER»ABS!!RHQB-RH0HB{IHS»SURF)) / RHOHBl IHS, SURF) ) ) 
60 CONTINUE 

IF (RELER.LT..001) I END= IEND+1 
WRITE(6, 1020) ITER.RELER 
MARK = 0 

WRITE ALL BLADE SURFACE VELOCITIES 

IF (SURVL.LE.O) RETURN 
WRI TE(6,1030) 

MBOT = MIN0(MMMM1,MB0C) 

I F ( M80T . L T. 2 ) GO TO 70 
WRI TE ( 6 , 1040 ) 

WRI TE16, 1050) I MV ( IM ) » WMB! IM , I ) , BE TAVC I M , 1 ) , WWCRM ( I M , l ) , WMB ( I M,2 ) , 
I BET AVI IM,2) , WWCRM ( IM,2), IM=2,MB0T) 

70 MBIT = MAXO ( 2 , MB I I ) 

IF(MBIT.GT.MMMMi) GO TO 80 
WRITEI6, 1060) 

WRITE! 6, 1050) l MV ( I M ) ,WM8( I M , 3 > , BETAV { I M , 3 ) , WWCRM l I M, 3 ) , WMB 1 1 M,4 ) , 
IBETAV! IM,4) , WWCRM! IM, 4), IM=MB I T , MMMM 1 ) 

80 WRITE! 6,1070) 

DO JO SURF= 1,4 
I MS S = IMS(SURF) 

IF l IMSS.EQ.O) GO TO 90 
WR I TE! 6, 1080 ) SURF 

WRI TE 16,1090) (MHt IHS, SURF ) ,WTB! IHS, SURF ) .BETAHI IHS, SURF ) .WWCRT 
1 ( IHS, SURF) , IHS=1, IMSS ) 

90 CONTINUE 
RETURN 

1000 FORMAT! 1H 1/// /40X , 34H VFLOC I T I ES AT INTERIOR MESH POINTS//) 

1010 FORMAT ! 1HL, 3HIM=, 13, 5 (24H VELOCITY ANGL E I DEG ) ) / 

1 !5X,5!G15.4, F9.2) ) ) 

1020 FORMAT! 14HLI TERATION NO. , I 3, 3X , 36HMAX IMUM RELATIVE CHANGE IN DENS1 
1TY = , G 1 1 . 4 ) 

1030 FORMAT ( 1H1////16X, 1H*»25X,49HSURFACE VELOCITIES BASED ON MERIDIONA 
1L COMPONENTS ,36X, IH* ) 

1040 FORMAT! 16X, IH* , 53X , 1 H* , 5 6X , 1H*/16X,1H*, 19X, 15HBLADE SURFACE 1,19X, 
11H*,20X, 15 H BLADE SURFACE 2,2 IX, 1H*/7X,1HM,8X,1H*,2(3X,8H VELOCITY, 
13X, 10H ANGLE! DEG) ,5X,5HW/ WCR, 19X, 1H*,3X) ) 

1050 FORMAT ! 1H ,G13.4,3H * , G 12 . 4, F9. 2, G l 5. 4, 1 7X , IH* , 3X , G12. 4 , F9. 2 , G1 5. 

14, 17X, IH* ) 

1060 FORMAT ! /// 16X, IH*, 19X , 15HBLADE SURFACE 3, 19X , IH* ,20X, 15HBLADE SURF 
1 ACE 4, 2 IX, 1H*/7X, 1HM, 8X, IH* , 2 l 3X , 8HVEL0C I TY , 3X , 10HANGLE I DEG ) , 5X , 
15HW/WCR,19X,1H*,3X)) 

1070 FORMAT! IH 1 // / / 3X , 49H SURF ACE VELOCITIES BASED ON TANGENTIAL COMPONE 
1NTS ) 

1080 FORMAT ( //22X, 15HBLADE SURFACE , I 1/7X, 1HM, 10X, 8HVEL0CITY,3X, 10HANG 
ILE(DEG) ,3X,5HW/WCR) 

1090 FORMAT I IH , 2G 1 3 . 4 , F9 . 2 , G 15 . 4 ) 

END 
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SUBROUTINE BLCD 


BLCD CALCULATES BLADE THETA COORDINATE AS A FUNCTION OF M 
COMMON SRW, ITER. I END , LER ( 2 ) , NER ( 1 ) 

COMMON /INP/ GAM , AR» T IP» RHOI P » WTFL * WTFL SP, OMEGA, ORF , BETAI ,BETAO, 
1NOBL.MBI ,MB0,MBI2,MB02,MM,NBBI,NBL,NRSP,MBDYF,MBDYL,ITF, ITL, 
2BLDAT, AANDK, ERSOR,STRFN, INTVL » SURVL » MAGF AC » 

3 MR ( 50) t RMSP ( 50 ) ,BESP(50) 

COMMON /CALCON/ MBI I , M800, MMM, MB 1 1 Ml , MB I I P 1 , MB00M1 , MB00P1 , MMMM1 , 
1HM1,HM2,HM3,HT,DTLR,DMLR,PITCH,CP,£XP0N,TWW.CPTIP,TGR0G,TBI,TB0, 
2LAMBDA.TWL, ITOR, ITMAX.NIP, IMS( 4 ) , BVl 4 ) , M V( 100) ,IV( 101) , 

3UBV( 100, A) ,RW6V( 100,4 ) , I TV ( 100, 6) , TV ( 100, 4 ) , DTDMV ( 1 00 , 4 ) , 

4BETAV1 100 ,4) ,MH( 100,4) , DTDMH UOO, 4 ) , BET AH { 100,4) ,RMH( 100,4) , 
5BEH{ 100,4) ,RM( 100) ,BE( 100) ,DBDM( 100) , SAL ( 100 ) , AAA ( 100 ) 

COMMON /GEOMIN/ CHORD l 4 ) , STGR ( 4 > ,MLE ( 4 ) , THLE ( 4 ) , RM I (4 ) , RMO ( 4 ) , 
1RI14),R0(4) , BET I ( 4 ) , BETO (4 ) , NSP II 4 ) , MSP I 50 , 4 ) ,THSP(50,4) 

COMMON / BLCDCM/ EM150.4) ,INITC4) 

ENTRY BL1(M, THETA, DTDM, INF) 

INTEGER BLDAT, AANDK, ERSOR, STRFN , SURVL , AATEMP , SURF , 

IF I RST, UPPER, SI, ST, SRW 

REAL K , KAK .LAMBDA, LMAX,MH,MLE, MR, MSP , MV.MVIM1 
REAL M,MMLE,MSPMM,MMMSP 
SURF= 1 
SI GN= 1. 

GO TO 10 

ENTRY BL2 (M, THETA, DTDM, INF) 

SURF= 2 
S I GN= -1. 

GO TO 10 

ENTRY BL3(M, THETA, DTDM, INF) 

SURF= 3 
SIGN= 1. 

GO TO 10 

ENTRY BL4(M, THETA, DTDM, INF) 

SURF= 4 
S I GN= -1. 

10 INF= 0 

NSP= NSPI (SURF) 

IF (INIT(SURF).EQ.13) GO TO 30 
IN I T ( SURF ) = 13 

INITIAL CALCULATION OF FIRST AND LAST SPLINE POINTS ON BLADE 

AA = BETI (SURFJ/57. 295779 
AA = SIN(AA) 

MSP ( 1 , SURF ) = RI (SURF)*(1.-SIGN*AA) 

BB = SORT ( 1 AA**2 ) 

THSP( 1 , SURF) = S IGN* BB*R I(SURF)/RMI(SURF) 

BET I l SURF ) = AA/BB/RMI (SURF) 

AA = BETO(SURF)/57. 295779 
AA = SIN(AA) 

MSP ( NSP , SURF ) = CHORD (SURF )— RO ( SURF ) * ( 1 . +S I GN* AA ) 

BB = SQRT ( l.-AA**2) 

THSP(NSP, SURF) = STGR ( SURF ) +S IGN*BB*RO ( SURF ) /RMO( SURF ) 
BETO(SURF) = AA/BB/RMOISURF) 

DO 20 I A= 1 , NSP 

MSP ( IA, SURF )= MSP ( IA, SURF) +MLE l SURF ) 

20 THSP ( I A , SURF ) = THSP ( I A , SURF ) +THLE ( SURF ) 

CALL SPLN22(MSP( 1, SURF), THSP l 1, SURF) , BETI (SURF) , BE TO (SURF) , NSP, 
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1AAA, EM{1, SURF) ) 

IF (BLOAT. LE.O) GO TO 30 
IF (SRW.E0.0) WRITEt 6,1000) 

SRW = 1 

WRITE(6, 1010) SURF 

WRITE (6,1020) (MSP( I A, SURF) ,THSP( IA,SURF) , AAA( IA) ,EM( IA.SURF) , 
1 1 A= 1 » MSP ) 

BLADE COORDINATE CALCULATION 
30 KK = 2 

IF (M.GT.MSP(ltSURF) ) GO TO 50 

AT LEADING EDGE RADIUS 

MMLE = M-MLE ( SURF ) 

IF (MMLE.LT.-DMLR) GO TO 90 
MMLE= AMAX1 (0. ,MMLE) 

THE T A= SQRT(MMLE*(2.*RI(SURF)-MMLE) )*SIGN 

IF ( THETA. EQ.O. ) GO TO 40 

RMM= R I ( SURF ) — MMLE 

DTDM= RMM/THETA/RMK SURF ) 

THETA= THETA/RMI (SURF )+THLE( SURF) 

RETURN 
40 I NF = 1 

DTDM = I . E 10*S IGN 
THET A= THLE ( SURF ) 

RETURN 

ALONG SPLINE CURVE 

50 IF (M.LE.MSP(KK, SURF) ) GO TO 60 
IF (KK.GE.NSP) GO TO TO 
KK = KK + 1 
GO 10 50 

60 S= MSP(KK,SURF)-MSP(KK— 1,SURF) 

EMKM 1= EM ( KK— 1 , SURF ) 

EMK= EM ( KK , SURF ) 

MSPMM= MS P ( KK , SURF ) -M 
M M M S P= M-MSP(KK-1,SURF) 

THK= THSP ( KK , SURF )/S 
THKMI= THSP(KK-1,SURF)/S 

THE TA= EMKMl*MSPMM**3/6./S + EMK*MMMSP**3/6. /S + ( THK-EMK*S/6. ) * 
1 MMMSP + (THKMI-EMKMl*S/6. )*MSPMM 
DTDM= -EMKMI*MSPMM**2/2. /S + EMK*MMMSP**2/2./S + THK-THKM1— ( EMK- 
1 EMKML ) *S/6. 

RETURN 

AT TRAILING EDGE RADIUS 

70 CMM= CHORDl SURF )+MLE ( SURFJ-M 
IF (CMM.LT.-DMLR) GO TO 90 
CMM= AM AX 1 ( 0 . , CMM ) 

THE TA= SORT (CMM*(2.*R0(SURF) —CMM ) )*SIGN 

IF ( THETA. EQ.O. ) GO TO 80 

RMM= RO ( SURF ) -CMM 

DTDM = -RMM/THETA/RMOISURF ) 

THETA = STGR(SURF)+THETA/RMO(SURF)+THLE(SURF) 

RETURN 
80 I NF= 1 
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DTDM = -1 . E 10*S IGN 
THETA= THLE(SURF)+STGR(SURF) 

RETURN 

ERROR RETURN 

90 WRITE (6,1030) LER ( 2 ) , M , SURF 
STOP 

1000 FORMAT ( 1H1 , 13X.27HBLADE DATA AT SPLINE POINTS) 

1010 FORMAT ( 1HL, 17X, 16HBLADE SURFACE, 14) 

1020 FORMAT (7X , 1HM, 10X, 5HTHETA, 10X, 10HDERIVATIVE, 5X, 10H2ND DERIV. / 

1 (4G15.5) ) 

1030 FORMAT (14HLBLCD CALL N0.,I3/33H M COORDINATE IS NOT WITHIN BLADE/ 
14H M =,G14.6, 10X.6HSURF =,G14.6) 

END 


FUNCTION I PF ( IM, IT) 

COMMON /CALCON/ MB 1 1 , MBOO, MMM, MBI I M 1 , MB I I P 1 , MB00M1 , MB OOP 1 , MMMM 1 , 
1HM 1,HM2.HM3,HT,DTLR,0 MLR, PITCH, CP, EXPON, TWW , CP T I P , TGROG , TBI ,TBO, 
2 LAMBDA, TWL, I TOR, ITMAX,NIP,IMS{4),BV(4),MV(100),IV{101), 

3UBV( 100,4) , RWB V ( 1 00 , 4 ) , I TV ( 1 00 , 6 ) , T V ( 100 , 4 ) , DTDMV 1 100 , 4 ) , 

48ETAVI 100 ,4) ,MH( 100,4 ) , DTDMH ( 100,4) , BET AH( 100,4) ,RMH( 100,4) , 
5BEH( 100,4) ,RM( 100) ,BE ( 100) ,DBDM( 100) , SAL ( 100) , AAA (100) 

IPF= I V ( I M ) + I T-2 

IF ( IT.LE. I TV( IM,4) ) RETURN 

I F ( IT.LE. ITV( IM,2).AND.IM.LT.MBII) RETURN 
I F ( IM.GE.MI NO (MBIIjMBOOPl) ) GO TO 10 
IF ( ITV( IM,1) .EQ.-10000) RETURN 
I PF= IPF-ITV(IM f l) + ITV(IM,2)+1 
RETURN 

10 I F ( IM.GT.MAXO( MB 1 1 Ml, MBOO) ) GO TO 20 
I F ( IM.GT.MROO) RETURN 
I PF= IPF— ITV( IM,3)+ITV( IM,4)+1 
I F ( IT.LT. ITV( IM, 1) ) RETURN 
I PF= IPF-ITV(IM f l)+ITV(IM,2)+l 
RETURN 

20 IF ( I T V ( IM,3).EQ.-10000) RETURN 
I PF= IPF-ITVl IM,3)+ITV( IM,4)+l 
RETURN 
END 


Subroutine BDVINT 

BDVINT calculates interpolated values of the stream function along either a vertical 
or a horizontal boundary for the fine mesh. The interpolation is based on a cubic spline 
curve (ref. 8) using the stream -function values at the original coarse-mesh points. 

The input arguments for BDVINT are as follows: 

BVIN input array from main program (see DESCRIPTION OF INPUT AND 

OUTPUT) 
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input array from main program (see DESCRIPTION OF INPUT AND 
OUTPUT) 

NSP given number of stream-function values on the coarse mesh 

MVTH array of either m(MV) or 0(TH) coordinates for the fine mesh 

KBDRY input variable (see DESCRIPTION OF INPUT AND OUTPUT) 

DLR tolerance, either DMLR or DTLR 

MMMITX total number of fine mesh lines , either MMM or ITMAX 

The internal variables for BDVINT are as follows: 

BVINT array of values of m or 6 where interpolated stream-function values 

are desired 

IM index of mesh line 

IMl value of IM for first point in region 

IM2 IMl + 1 

IMT value of IM for last point in region 

NS PENT number of interpolated values of stream function 

RWBV array of interpolated values of either 9u/3m or du/39 

UBV array of interpolated stream-function values 

SUBROUTINE BDV I NT ( BV I N , U BV I N , NSP , MVTH, KBDR Y , DLR , MMMI TX ) 

COMMON /CALCON/ MBII , M800, MMM, MB 1 1 Ml , MB I I P 1 , MB00M1 , MBOOP l ,MMMM1 , 
1HM1 , HM2.HM3.HT »DTLR» DMLR»PITCH»CP»EXPQN,TWW»CPTIP» TGROG, TBI ,TBO, 
2LAMBDA,TWL» IT OR, ITMAX, NIP, IMS!4),BV(4),MV{ 100) , IV 1 101 ) , 

3UBVI 100,4 ) .RWBV ( 100, 4) , I TV { 100,6) , TV! 100 , 4 ) , OTOM V ( 100, 4 ) , 

4BETAV! 100,4) ,MH( 100, 4 ) , DTDMH ( 1 00, 4 ) , BET AH { 100, 4 ) , RMH ( 100 , 4 ) , 

5BEHI 100,4) ,RM{ 100) ,BEI 100) ,DBDM( 100) ,SAL ( 100) , AAAI 100) 

DIMENSION BVINT I 100) ,BVIN{ 100 ), UBV I N I 100 ), MVTH 1 100 ) 

REAL MVTH 

DO 10 IM=2, MMMITX 

IF (MVTH! IM) .GT.BVIN I 1)+DLR) GO TO 20 
10 CONTINUE 
IM = MMMITX 
20 IMl = IM-1 

BVINT! IMl ) = BVIN(l) 

I M2 = IM 

DO 30 I M= I M2, MMMITX 

IF I MVTH! IM).GT.BVIN(NSP)-DLR) GO TO 40 
30 BVINT! IM) = MVTH ( I M) 

IM = MMMITX 
40 I MT= I M 

BVINT! IMT) = BVIN1NSP) 

NSP I NT= IMT-IM1+1 

CALL SPLINT! BV IN, UBV IN, NSP, BVINT! I M 1 ) , NSPI NT, UBV ( I Ml , KBDRY ) , RWBV ( I 
IMl, KBDRY) ) 

RETURN 

END 
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Program Listing of Remaining Subroutines 


The remaining subroutines are described in reference 3. The description applies 
even though the subroutines have been revised. 


SUBROUTINE MHORIZIMV, I TV ,BL, MB I ,MBO, I TO , HT, DTLR , KODE , J , MH, DTDMH, 
IMRTS) 

C 

C MHORIZ CALCULATES M COORDINATES OF INTERSECTIONS OF ALL HORIZONTAL 
C MESH LINES WITH A BLADE SURFACE 
C KODE = 0 FOR UPPER BLADE SURFACE 
C KODE = 1 FOR LOWER BLADE SURFACE 
C 

COMMON SRW, ITER, I END , LER !2 ) , NER !1) 

DIMENSION MV ( 100) , IT V I 100 ) , MH ! 100 ) .DTDMH! 100) 

INTEGER BLDAT , AANDK* ERSOR, STRFN , SURVL , AATEMP , SURF , 

IF I RST, UPPER, SI, ST, SRW 

REAL K,KAK, L AMBDA, LM AX , MH, ML E, MR , MSP , MV, MV I M 1 
REAL MV I M 
EXTERNAL BL 

IF IM8I.GE.MB0) RETURN 
IM= MBI 
10 I TI ND= 0 

20 IF ( ITV( IM+l) — ITVI IMJ-IT IND) 30,40,50 
30 J= J+l 

T I = FLOAT ( ITVI I M+ 1 )- I TO- IT I ND+KODE ) *HT 
I T I ND= ITIND-1 
MV I M = MV I IM) 

IF (MRTS.EO.l) MV I M = MV IM + I MV I IM+ l )-M V IM ) / 1 000. 

CALL ROOT !MVIM,MV!IM+1) ,TI,BL, DTLR, MHIJ), DTDMH! J) ) 

GO TO 20 
40 I M= IM+l 
MRTS = 0 

IF (IM.EQ.MBO) RETURN 
GO TO 10 
50 J= J+l 

T I = FLOAT! ITVI I M )- ITO+ I T IND + KODE ) *HT 
I TI ND= ITIND+1 
MVIM = MV! IM) 

IF IMRTS. EO.l) MVIM = MV IM+ ( MV ! IM+ 1 )-MV I M ) / 1 000. 

CALL ROOTIMVIM , MV! I M+ 1 ) , T I , BL , DTLR , MH I J ) , DTDMH ! J ) ) 

GO TO 20 
END 
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SUBROUTINE DENSTY ( RHOW , RHO , VEL , TWLMR , CPT IP , EXPON ,RHOI P , GAM, AR , TI P ) 

DENSTY CALCULATES DENSITY AND VELOCITY FROM THE WEIGHT FLOW PARAMETER 
DENSITY TIMES VELOCITY 

COMMON SRW»ITER, IEND, LERI2) »NER( l) 

VEL = RHOW/RHO 
IF (VEL.NE.O.) GO TO 10 
RHO = RHO IP 
RETURN 

10 TTIP = l.-(VEL**2+TWLMR)/CPTIP 
IF(TTIP.LT .0. ) GO TO 30 
TEMP = TT I P** ( EXPON— 1 . ) 

RHOT = RHOIP*TEMP*TT IP 

RHOWP = -VEL **2/ GAM* RHO IP/AR*TEMP/TI P+RHOT 
I F ( RH0WP.L6.0. ) GO TO 30 
VELNEW = VEL+ ( RHOW— RHOT * VEL ) /RHOWP 
IF ( ABS ( VEL NEW— VEL )/VELNEW.LT.. 0001 ) GO TO 20 
VEL = VELNEW 
GO TO 10 
20 VEL = VELNEW 
RHO = RHOW/VEL 
RETURN 

30 TGROG = 2.*GAM*AR/ (GAM+1.) 

VEL = SQRT(TGROG*TIP*( l.-TWLMR/CPTIP) ) 

RHO = RH0IP*(1.-(VEL**2+TWLMR)/CPTIP)**EXP0N 
RWMORW = RHOW/RHO/VEL 
NER(l) = NER ( 1 ) + l 

WRITE (6,1000) LER(I ) »NER( 1) , RWMORW 
IF (NER( D.EQ.50) STOP 
RETURN 

1000 FORMAT ( 1 6HLDENSTY CALL N0..I3/9H NER ( 1 ) =,I3/10H RHO*W IS ,F7.4, 
134H TIMES THE MAXIMUM VALUE FOR RHO*W) 

END 


SUBROUTINE ROOT ( A, B, Y , FUNCT , TOLER Y, X , DFX ) 

C 

C ROOT FINDS A ROOT FOR (FUNCT MINUS Y) IN THE INTERVAL ( A, B ) 
C 

COMMON SRW, ITER, I END , L ER ( 2 ) , NER ( I) 

INTEGER SRW 

IF (SRW.EQ.21) WRITE16.1000) A , B , Y , TOLER Y 
TOLERX= ( B— A ) /1000. 

AB2= (A+B)/2. 

1= 0 
X= A 

10 CALL FUNCT(X,FX,DFX, INF) 

IF (SRW.E0.21) WR I TE ( 6, 1010 ) I , X, FX, DFX , INF 
IF (ABS(Y-FX).LT.TOLERY) RETURN 
IF (I.GE.1000) GO TO 30 
1= I+I 

IF ( INF.NE.O .OR. DFX.EQ.O.) GO TO 20 
X= ( Y-FX ) / DFX+X 

IF (X.GE.A .AND. X.LE.B) GO TO 10 
X = A+TOLERX*FLOAT( I ) 
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IF ( I .60. 1 ) X = B 
GO TO 10 

20 IF IX.LT.A82) X=X+TOLERX 
IF (X.GE.A82) X=X— TOL ERX 
GO TO 10 

30 WRI IE I 6 , 1 020 ) LER(2),A,B,Y 
STOP 

1000 FORMAT ( 32H1 INPUT ARGOMENTS FOR ROOT — A =G13.5,3X» 3HB =,G13.5, 
13X.3HY =»G13.5«3X*8HT0LEKY =,G13.5/17H ITER. 00. X,17X, 

22HFX,15X,3HDFX,10X,3HINF ) 

1010 FORMAT (5X, I 3 , G 16. 5, 2G 18.5 , 16) 

1020 FORMAT (14HLR00T CALL N0.,I3/47H ROOT HAS FAILED TO CONVERGE IN 10 
100 I TERATI0NS/4H A = , G14 .6 , 1 OX, 3HB = , G14. 6 , 1 OX , 3HY =,G14.6) 

END 


SUBROUTINE SPLINE (X, Y,N, SLOPE, EM ) 

SPLINE CALCULATES FIRST AND SECOND DERIVATIVES AT SPLINE POINTS 
END CONDITION - SECOND DERIVATIVES ARE THE SAME AT END POINT AND 
ADJACENT POINT 

COMMON Q/BOX/S ( 100 ), A ( 100) ,B( 100) ,C( 100 ) , F I 100 ) , W ( 100 ) , SB ( 1 00 ) , 
1GC200) 

DIMENSION X(N) *Y(N),EM(N),SLQPE(N) 

INTEGER Q 
DO 10 1=2, N 
10 S ( I )=X( I ) -X ( I- 1 ) 

N0= N— 1 

I F ( NO. L T. 2 ) GO TO 30 
DO 20 1=2, NO 
A(I)=S(I)/6. 

B< I ) = t SI I )+S ( 1 + 1 ) )/3. 

C(I)=S I 1+1 )/6. 

20 F< I ) = { YII + 1 )-Y( I ) )/S( 1 + 1 )-{ Y( I )— Y C 1-1) )/S< I ) 

30 AIN) = -.5 
B ( 1 ) =1 . 

B ( N ) = 1 • 

Cl 1) = -.5 
F ( 1 ) =0 . 

F ( N ) =0 . 

WI 1)=B( 1) 

SB ( 1 ) =C ( 1 ) / W ( 1 ) 

G 1 1 ) =0 . 

DO 40 1=2, N 

W ( I )=B( I )— AII)*SB{ I-l ) 

SB ( I )=C 1 1 ) /W ( I ) 

40 G( I )=(F(I )-A(I)*G( I-l ) )/W(I) 

EMC N )=G ( N ) 

DO 50 1=2, N 
K=N+ 1— I 
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50 EM(K)=G(K)-SB(K)*£M(K+l) 

SLOPE l 1)=— S(2)/6.*(2.*EM(1)+EM{2) )+(Y(2)-Y(l))/S(2) 

DO 60 I =2 , N 

60 SLOPE ( I) = S(I )/6.*C2.*EMC D+EHC 1-1) ) + (Y(I)-Yt l-l))/sm 

IF (Q.EQ.13) WRITE (6,1000) N, ( X ( I) , Y< I ) , SLOPE ( I) , EM (I) , 1=1 , N) 
RETURN 

1000 FORMAT (2X,15HN0. OF POINTS =, I 3/10X , 1HX, 19X , 1HY, 1 9X, 5HSL0PE , 1 5X , 
12F)EM/(4F20.8) ) 

END 


SUBROUTINE SPLN 22 ( X , Y , Y IP , YNP , N, SLOPE , EM > 

SPLN 22 CALCULATES FIRST AND SECOND DERIVATIVES AT SPLINE POINTS 
END CONDITION - DERIVATIVES SPECIFIED AT END POINTS 

COMMON Q/BOX/S ( 100 ). A ( 100 ) , B I 100 ) , C ( 100 ) ,F ( 100 ) , W { 100 ) , SB t 100 ) , 
1G ( 200) 

DIMENSION X(N) ,Y(N),EM(N),SLOPE(N) 

INTEGER 0 
DO 10 1 = 2 , N 
10 S( 1 )=X( I ) — XI I — 1 ) 

NO=N — 1 

IFIN 0 .LT. 2 ) GU TO 30 
DO 20 1 = 2 , NO 
A ( I ) =S ( I ) / 6 . 

B( I ) = (S( I )+S< 1 + 1 ) )/ 3 . 

C C I )=S( H- 1 J/ 6 . 

20 F(I) = (Y(I + l)-Y(I))/S(I + l)-(Y(I)-Y(I-l))/Sm 
30 A ( N ) = S ( N ) / 6 . 

B (1) =S ( 2 ) / 3 . 

BIN) = S(N)/ 3 . 

C(l)=S( 2 )/ 6 . 

F( 1 ) = (Y( 2 ) — Y(l) )/S( 2 )-YlP 
FIN) = YNP-IY(N)— Y(N— 1 ) ) /SIN) 

W I 1 ) =B I 1 ) 

SB ( l ) =C I 1 )/Wll) 

G ( 1 )=F( 1 ) /WI 1 ) 

DO 40 1 = 2 , N 

W( I )=B l I ) —At I ) *SB ( 1-1 ) 

SB ( l ) =C ( I ) /W ( I ) 

40 G( I ) = (F(I )-A(I)*GI 1 - 1 ) )/Wl I) 

EM ( N ) =G ( N ) 

DO 50 1 = 2 , N 
K=N+ 1 -I 

50 EM(K)=G(K)-SB(K)*EM(K+ 1 ) 

SLOPE ( 1 )=-SI2)/6.*I2.*EM(l )+EMI2) ) + 1 Y( 2 )-Y 1 1 ) ) /S (2 ) 

DO 60 1 = 2 , N 

60 SLOPE! I ) = S I I )/ 6 . * 12 .* EMI I)+EMI 1 - 1 ) ) + {Y( I ) — Y I I — 1 ) ) / S I I ) 

IF (Q.EQ. 18 ) WRITE ( 6 , 1000 ) N, I X ( I ) , Y 1 1 ) , SLOPE (I) , EM (I) , 1 = 1 , N ) 
RETURN 

1000 FORMAT ( 2 X. 15 HN 0 . OF POINTS =, I 3 / 10 X , 1 HX, 19 X , IHY , 1 9 X , 5 HSL 0 PE , 15 X , 
12 HEM/( 4 F 20 . 8 ) ) 

END 
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SUBROUTINE SPLINT ( X , Y ,N ,Z , MAX, YI NT , DYDX ) 

C 

C SPLINT CALCULATES INTERPOLATED POINTS AND DERIVATIVES FOR 
C A SPLINE CURVE 

C END CONDITION - SECOND DERIVATIVES ARE THE SAME AT END POINT AND 
C ADJACENT POINT 
C 

COMMON Q/BOX/St 100), A( 100) , B { 100 ) • C ( 1 00 ) ,F ( 100 ) , W ( 1 00) , SB ( 100 ) , 
lG(lOO) ,EM(100) 

DIMENSION X(N),Y(N) , Z < MAX ). Y INT { MAX ), DYDX ( MAX ) 

INTEGER Q 

IF C MAX . LE • 0 ) RETURN 
III = Q 
DO 10 1=2, N 

io sm=xm-xu-i) 

NO=N— 1 

IFINO.LT. 2) GO TO 30 

DO 20 1=2, NO 

All )=S< I )/6.0 

Bl I ) = (S(I )+Sl 1 + 1) )/3.0 

ci i ) =s i i + n/6.0 

20 Fll ) = l Yl 1 + 1 )-Yl 1 ) J/S l 1 + 1 )-l YID-Yl 1-1 ) )/Sl I ) 

30 AIN) = -.5 
Bl 1 )=1.0 
B t N ) =1 . 0 

cm = -.5 

FI 1 )=0.0 
FIN) =0 . 0 
Wll )=B( 1) 

SB 1 1 ) = C l 1 )/W(l) 

G l 1 ) =0 . 0 
DO 40 I =2 , N 

Wll )=Bl I )— Al I )*SB( 1-1 ) 

SBl I )=C( I )/WlI) 

40 Gl I ) = l F l I > -A l I )*Gl 1-1 ) )/Wl I ) 

EM t N) = G ( N ) 

DO 50 1=2, N 
K=N+1— I 

50 EM(K)=GIK)-SBIK)*EMIK+1) 

DO 140 1=1, MAX 
K=2 

IF IZl I ) -X ( 1) ) 70,60,90 
60 YINTII ) =Y 11) 

GO TO 130 

70 I F l Z ( I ) . GE . 1 1 . 1*X ( 1 ) — . l*Xl 2 ) )) GO TO 120 
WRITE 16,1000) Zll) 

0 = 16 
GO TO 120 
80 K=N 

IFIZID.LE. ll.l*X(N)-.l*X(N-l) ) ) GO TO 120 
WRITE 16,1000) Zll) 

Q = 16 
GO TO 120 

90 IF l Z l I ) -X l K ) ) 120,100,110 



100 Y I N T { I } = Y ( K ) 

GO TO 130 
110 K=K + 1 

IF(K-N) 90,90,80 

120 YINT(I) = EM(K-1)*(X(K)-Z( I) > **3/6 . / S ( K J +EM ( K ) * ( Z { I J-XlK-l) )**3/6. 
1/S (K)+ ( Y(K) /S(K)-EM( K ) *S {K ) / 6. } * ( Z ( I ) -X ( K~ 1 ) ) + I Y < K- l ) / S I K) -EM { K- 1 ) 
2*S(K)/6. > * ( X C K » — Z ( I) ) 

130 DYDXU ) =— EM ( K— 1 ) * ( X { K ) — Z ( I ) ) **2/2 . 0/S ( K ) +EM ( K ) * ( X ( K- 1 )-Z (I) )**2/2. 

10/S (K) + ( YIK)-Y(K-l) )/S(K >-(EM(K)-EM(K-l) )*S(K)/6.0 
140 CONTINUE 

MX A = MAX 0 ( N » MAX ) 

IF ( U.EG. 16 ) WRITE (6, 1010) N,MAX, I X { I ) , Y ( I ) , Z ( I ) , YINTCI) , DYDX (I) , 
11=1 ,MXA) 

Q = III 
RETURN 

1000 FORMAT (54H SPLINT USED FOR EXTRAPOLATION. EXTRAPOLATED VALUE = , 
1G14.6) 

1010 FORMAT ( 2X , 2 1HN0. OF POINTS GIVEN =,I3,30H, NO. OF INTERPOLATED PO 
1 1 NTS =, I3/10X, 1HX, 19X, 1HY, 16X, 11HX-INTERP0L. ,9X,11HY- INTERPOL. , 
28X, L4HDYDX- INTERPOL. / ( 5E20.8 ) ) 

END 


Lems Research Center, 

National Aeronautics and Space Administration, 
Cleveland, Ohio, December 12, 1968, 
126-15-02-31-22. 
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